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Abstract 

This  thesis  develops  a  method  for  estimating  the  normal  mode  decomposition  of 
broadband  signals  and  uses  it  to  analyze  data  from  the  Acoustic  Thermometry  of 
Ocean  Climate  (ATOC)  experiment.  Normal  modes  are  the  eigenfunctions  of  the 
ocean  waveguide,  derived  from  the  frequency-domain  wave  equation.  They  are  useful 
in  underwater  acoustics,  particularly  matched  field  processing  and  tomography,  be¬ 
cause  the  lowest  modes  provide  an  efficient  description  of  the  most  energetic  arrivals 
at  long  ranges.  Extracting  source  or  environmental  information  from  the  mode  sig¬ 
nals  depends  on  understanding  the  effects  of  internal  waves  on  coherence  and  the 
validity  of  adiabatic  approximation.  While  much  theoretical  research  has  been  done 
on  long-range  propagation  of  modes  in  deep  water,  there  have  been  few  opportunities 
to  compare  theoretical  predictions  with  experimental  measurements. 

The  first  contribution  of  this  thesis  is  a  short-time  Fourier  framework  for  es¬ 
timating  broadband  signals  propagating  in  the  lowest  modes  of  the  ocean  wave¬ 
guide.  Since  previous  research  has  focused  primarily  on  narrowband  sources,  this 
work  concentrates  on  broadband  processing  issues.  Specifically,  it  addresses  the  fun¬ 
damental  issue  of  frequency  resolution  required  for  mode  estimation,  analyzes  the 
performance  characteristics  of  two  modal  beamforming  algorithms  and  explores  the 
time/frequency  tradeoffs  inherent  in  STFT  mode  processing. 

1  he  second  contribution  of  this  research  is  a  detailed  analysis  of  the  low-mode 
arrivals  at  megameter  ranges  using  five  months  of  data  from  the  ATOC  vertical  line 
array  at  Hawaii  (3515  km  range).  Short-time  Fourier  processing  of  these  receptions 
revealed  that  each  low  mode  contains  a  series  of  arrivals,  rather  than  the  single  dis¬ 
persive  arrival  that  would  characterize  adiabatic  propagation.  Average  coherence 
times  of  the  mode  signals  are  on  the  order  of  6-8  minutes.  The  multipath  structure 
changes  significantly  between  receptions  at  4-hour  intervals,  indicating  that  stochas- 


2 


tic  methods  are  required  for  mode  tomography  at  megameter  ranges.  A  statistical 
analysis  found  that  modes  do  retain  travel-time  information  at  megameter  ranges, 
e.g the  centroids  show  the  expected  dispersion  characteristics  of  a  deep  water  chan¬ 
nel.  The  centroids  show  statistically  significant  trends  in  mode  arrival  time  over  the 
period  of  the  experiment. 

Thesis  Supervisor:  Arthur  B.  Baggeroer 
Ford  Professor  of  Engineering 

Secretary  of  the  Navy/Chief  of  Naval  Operations  Chair  for  Ocean  Science 
Massachusetts  Institute  of  Technology 

Thesis  Supervisor:  James  C.  Preisig 

Assistant  Scientist,  Woods  Hole  Oceanographic  Institution 
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Chapter  1 


Introduction 


Normal  modes  provide  a  convenient  description  of  low-frequency  sound  in  the  deep 
ocean.  Their  strong  connection  to  the  propagation  environment  makes  them  useful 
in  a  variety  of  applications,  including  source  localization  and  acoustic  tomography. 
Currently,  there  is  much  interest  in  using  modes  to  analyze  broadband  receptions 
at  megameter  ranges  for  the  purpose  of  studying  ocean  variability  on  basin-scales. 
At  these  ranges,  the  effects  of  internal  waves  on  mode  coherence  are  not  known. 
This  thesis  develops  a  signal  processing  framework  for  estimating  modal  time  series 
and  uses  it  for  analyzing  data  from  the  Acoustic  Thermometry  of  Ocean  Climate 
experiment.  From  a  signal  processing  perspective,  the  key  issue  to  consider  is  the 
broadband  nature  of  the  signals;  specifically,  any  approach  must  accommodate  vari¬ 
ations  in  the  mode  characteristics  across  the  bandwidth  of  the  source.  In  examining 
the  data,  the  focus  is  on  understanding  the  fluctuations  of  mode  arrivals  and  char¬ 
acterizing  the  complicated  multipath  structure. 

The  rest  of  this  chapter  introduces  the  research  questions  addressed  by  this  thesis. 
As  a  starting  point,  the  following  section  motivates  the  use  of  the  modal  description  in 
the  context  of  long-range  acoustics  and  discusses  open  questions  about  megameter 
propagation.  The  second  section  highlights  some  of  the  signal  processing  issues 
surrounding  the  broadband  mode  estimation  problem.  Section  1.3  describes  the 
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Figure  1-1:  Deep  ocean  waveguide.  The  left  panel  shows  a  typical  deep  water  sound 
speed  profile.  The  right  panel  illustrates  how  refractive  effects  permit  propagation 
over  extremely  long  ranges. 

opportunities  presented  by  the  recent  ATOC  experiment.  Finally,  Section  1.4  states 
the  specific  research  objectives  and  outlines  the  remainder  of  the  thesis. 


1.1  Use  of  Modes  in  Long-Range  Acoustics 

The  deep  ocean  is  an  efficient  channel  because  it  traps  low-frequency  acoustic  signals, 
enabling  them  to  be  detected  thousands  of  kilometers  from  their  source.  Figure  1-1 
shows  how  the  refractive  effects  of  the  underwater  waveguide  allows  propagation  to 
such  long  ranges.  Acoustically,  the  deep  ocean  is  characterized  by  a  sound  speed 
profile  with  a  minimum,  between  800  and  1200  meters  depth,  known  as  the  sound 
channel  axis.  Sound  waves  bend  towards  regions  of  lower  velocity,  thus  the  minimum 
creates  a  duct.  As  the  figure  depicts,  a  sound  wave  leaving  the  source  on  a  downward 
trajectory  bends  back  towards  the  axis.  Once  it  passes  through  the  minimum  on  an 
upward  path,  it  bends  away  from  the  surface.  Purely  refracted  paths,  such  as  the 
one  shown,  do  not  scatter  energy  at  boundary  interactions.  Since  absorption  losses 
for  low  frequencies  (on  the  order  of  100  Hz)  are  minimal,  the  signals  can  propagate 
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Figure  1-2:  Broadband  reception  on  40-element  VLA  located  3515  km  from  source 


over  extremely  long  distances  in  the  channel. 

Fig.  1-2  illustrates  some  general  features  of  pulse  propagation  to  megameter 
ranges.  The  plot  shows  the  pulse-compressed  and  sequence-averaged  time  series 
recorded  by  a  40-element  vertical  line  array  (VLA)  located  3515  km  from  a  broad¬ 
band  source.  Inter-element  spacing  is  35  m,  and  the  array  is  approximately  centered 
on  the  sound  channel  axis.  This  figure  reveals  an  important  characteristic  of  propa¬ 
gation  in  the  underwater  sound  channel,  namely  the  time-spread  of  the  arrivals  due 
to  the  fact  that  they  take  many  different  paths  between  source  and  receiver.  The 
early  arrivals  traverse  deep-diving  ray  paths  that  are  associated  with  higher  group  ve¬ 
locities  because  they  sample  the  water  away  from  the  sound  speed  minimum.  Signals 
that  propagate  almost  horizontally,  along  the  sound  channel  axis,  arrive  last  and  are 
often  more  energetic.  In  general,  the  multipath  arrival  structure  can  be  represented 
in  terms  of  the  vertical  eigenfunctions,  or  normal  modes,  of  the  underwater  wave¬ 
guide.  Modal  dispersion  accounts  for  the  time-spread  of  the  signal  at  long  ranges.  In 
deep  water,  the  high  modes  travel  faster  than  the  low  modes,  thus  the  high  modes 
are  associated  with  the  early-arriving  energy  in  Fig.  1-2.  The  planewave-type  arrivals 
visible  in  the  early  parts  of  the  reception  (up  until  «  2373  seconds)  are  the  result 
of  constructive  interference  of  groups  of  these  higher  order  modes.  Planewaves  are 
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not  evident  in  the  last  2.5  seconds  of  the  reception,  associated  with  the  low  mode 
arrivals. 

As  solutions  to  the  frequency-dependent  wave  equation,  the  modes  provide  many 
useful  insights  about  sound  propagation.  Each  mode  essentially  samples  a  different 
part  of  the  water  column:  the  low  modes  are  concentrated  around  the  axis,  whereas 
the  higher  modes  have  greater  vertical  extents.  Since  the  modes  depend  strongly 
on  the  environment,  they  can  be  used  as  observables  for  matched  field  processing 
or  tomography  applications.  The  key  to  matched  field  or  tomographic  inversions 
is  the  ability  to  associate  an  arrival  with  a  particular  path  or  section  of  the  water 
column.  In  range-invariant  environments,  this  problem  is  trivial  because  the  modes 
propagate  independently  without  exchanging  energy,  i.e.,  an  arrival  in  mode  1  is 
known  to  have  traversed  the  entire  path  in  mode  1.  For  a  realistic  ocean  environment, 
however,  inhomogeneities  cause  coupling  of  energy  among  the  modes,  which  makes 
the  problem  much  more  difficult.  Understanding  the  mechanisms  and  effects  of  mode 
coupling  is  crucial  to  using  these  signals  in  any  type  of  application. 

At  long  ranges,  internal  waves  are  thought  to  be  the  primary  source  of  coupling. 
Vertical  displacements  of  water  associated  with  internal  waves  cause  fluctuations  of 
the  temperature,  thus  changes  in  the  sound  speed,  at  a  fixed  depth.  The  horizontal 
variability  of  these  sound  speed  fluctuations,  in  turn,  can  cause  an  exchange  of  en¬ 
ergy  among  the  modes  as  they  propagate.  The  effects  of  these  fluctuations  on  the 
planewave-type  arrivals  are  fairly  well-understood  -  these  arrivals  are  amenable  to 
analysis  via  geometrical  optics  approximation.  Significantly  less  is  known  about  the 
axial  mode  arrivals  since  there  is  no  comparable  theory.  The  late-arriving  modes 
tend  to  describe  the  most  energetic,  trapped  signals  and  thus  are  useful  in  detect¬ 
ing/estimating  weak  sources  at  long  range.  To  develop  an  understanding  of  how  these 
low  modes  propagate  through  internal  wave  fields  and  to  test  some  of  the  limited 
theoretical  results,  it  is  necessary  to  study  them  experimentally. 
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1.2  Broadband  Mode  Estimation  at  Megameter 
Ranges 

Measuring  the  mode  arrival  structure  at  long  ranges  presents  an  interesting  signal 
processing  problem.  Unlike  the  planewave  arrivals,  axial  mode  arrivals  generally 
overlap  in  time,  and  must  be  estimated  via  spatial  processing.  Since  the  modes  are 
an  orthonormal  basis  in  depth,  in  principle  they  can  be  separated  using  vertical  line 
arrays  spanning  the  entire  water  column.  In  practice,  the  degree  of  orthogonality 
of  the  modeshapes,  as  sampled  by  a  practical  array,  determines  how  well  the  modes 
can  be  resolved. 

A  key  issue  in  this  thesis  is  the  use  of  broadband  signals.  The  modes  are  inher¬ 
ently  frequency-dependent,  since  they  are  derived  from  the  frequency-domain  wave 
equation.  Previous  work  on  mode  estimation  has  primarily  focused  on  situation- 
s  where  a  narrowband  approximation  is  valid,  i.e.,  either  the  source  is  CW  or  the 
mode  functions  are  approximately  constant  across  the  band  of  the  source.  A  few 
researchers  have  implemented  broadband  mode  processors  using  an  FFT  for  the  fre¬ 
quency  decomposition,  but  they  have  not  discussed  the  frequency  resolution  required 
for  this  approach.  What  is  needed  is  a  general  framework  for  broadband  mode  esti¬ 
mation  that  will  allow  a  careful  analysis  of  performance  in  terms  of  mode  resolution 
and  time/frequency  resolution. 

A  recent  experiment  provides  an  opportunity  to  develop  methods  of  mode  pro¬ 
cessing  and  to  apply  them  to  studying  the  coherence  of  mode  arrivals  at  megameter 
ranges. 
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Figure  1-3:  ATOC  source  and  receivers 


1.3  Acoustic  Thermometry  of  Ocean  Climate 
Experiment 

The  purpose  of  the  Acoustic  Thermometry  of  Ocean  Climate  (ATOC)  experiment 
is  to  study  long-range  propagation  of  sound  and  to  investigate  acoustic  methods 
for  monitoring  ocean  climate  variability.  The  intent  is  to  demonstrate  that  travel¬ 
time  tomography  can  be  used  to  measure  ocean  temperature  over  ranges  of  3,000  to 
10,000  km.  The  ATOC  network  consists  of  a  broadband  source  off  the  California 
coast,  two  vertical  line  arrays,  and  a  number  of  bottom-mounted  horizontal  arrays. 
This  thesis  focuses  on  analyzing  data  from  the  two  VLA’s,  which  were  designed  to 
spatially  resolve  the  lowest  10  modes  at  each  location.  Figure  1-3  shows  the  location 
of  source  and  receivers  considered  in  this  thesis.  The  path  lengths  to  the  receivers 
at  Hawaii  and  Kiritimati  (Christmas  Island)  are  3515  km  and  5171  km,  respectively. 


19 


These  arrays  were  deployed  in  November  1995  and  recovered  in  September  of  1996. 
The  bottom-mounted  source  (934  m)  on  Pioneer  Seamount  transmitted  pulses  at 
a  center  frequency  of  75  Hz.  Each  transmission  consists  of  40  periods  of  a  pseudo¬ 
random  sequence,  phase  modulated  onto  a  75  Hz  carrier.  The  receiver  averages  every 
4-periods  internally.  Over  the  duration  of  the  experiment,  transmissions  are  sent 
every  4  hours  during  periods  established  by  the  ATOC  Marine  Mammal  Research 
Program. 

ATOC  presents  the  first  opportunity  to  study  mode  arrivals  at  megameter  ranges. 
The  next  section  outlines  the  objectives  of  this  thesis. 

1.4  Thesis  Objectives 

The  first  objective  of  this  research  is  to  define  a  framework  for  broadband  mode 
estimation.  Since  the  modeshapes  are  frequency-dependent  and  the  mode  spectral 
coefficients  are  time-dependent,  mode  estimation  involves  a  combination  of  temporal 
and  spatial  filtering.  Most  previous  work  has  focused  primarily  on  the  narrowband 
mode  estimation  problem,  and  has  not  addressed  issues  unique  to  broadband  signals. 

The  second  objective  is  to  analyze  the  low-order  mode  arrivals  in  the  ATOC  data. 
This  is  really  the  first  opportunity  of  its  kind.  Specifically,  this  research  hopes  to 
characterize  the  complicated  mode  arrival  structure  and  explore  the  effects  of  internal 
waves  on  mode  coherence. 

The  organization  of  the  thesis  is  as  follows.  Chapter  2  reviews  background  about 
normal  mode  representations  and  motivates  several  specific  questions  about  long- 
range  mode  propagation.  It  clearly  formulates  the  broadband  mode  estimation 
problem  and  proposes  and  approach  for  exploring  the  scientific/signal  processing 
research  topics  using  the  ATOC  data.  Following  that,  Chapter  3  presents  a  frame¬ 
work  for  broadband  mode  processing,  based  on  short-time  Fourier  techniques.  The 
fourth  chapter  presents  an  analysis  of  the  ATOC  data  set  for  the  Hawaii  array, 
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compares  these  results  to  simulations,  and  identifies  several  useful  statistics  of  the 
mode  arrivals.  Finally,  Chapter  5  summarizes  the  thesis  contributions  and  indicates 
directions  for  future  research. 
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Chapter  2 


Background 


As  indicated  in  Chapter  1,  normal  modes  are  of  interest  in  applications  such  as 
tomography  and  thermometry  because  the  lowest  modes  provide  a  convenient  de¬ 
scription  of  the  energetic  late  arrivals  at  megameter  ranges.  This  chapter  lays  the 
groundwork  for  the  rest  of  the  thesis  by  motivating  specific  questions  about  long- 
range  mode  propagation,  clearly  formulating  the  broadband  mode  estimation  prob¬ 
lem,  and  proposing  an  approach  for  exploring  these  research  topics  using  data  from 
the  ATOC  vertical  arrays.  The  material  is  divided  into  four  parts.  Section  2.1  re¬ 
views  the  salient  characteristics  of  the  normal  mode  representation  and  outlines  some 
of  the  open  questions  concerning  mode  propagation  through  range-dependent  and 
random  environments.  In  the  course  of  describing  relevant  features  of  the  modal 
basis  set,  this  section  also  introduces  a  range-dependent  ocean  environment,  which 
is  used  for  many  of  the  examples  in  later  chapters.  Given  this  mathematical  back¬ 
ground  and  experimental  motivation,  Section  2.2  poses  the  mode  estimation  problem 
for  vertical  arrays  and  highlights  important  design  considerations.  In  particular,  the 
discussion  emphasizes  the  broadband  character  of  the  problem  since  most  prior  work 
has  focused  on  using  modes  to  analyze  narrowband  signals.  The  third  section  reviews 
previous  work  on  mode  estimation  in  order  to  place  the  current  research  in  context. 
Finally,  Section  2.4  outlines  the  proposed  approach,  which  is  based  on  short-time 
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Fourier  analysis  techniques,  and  describes  the  areas  addressed  by  the  rest  of  the 
thesis. 


2.1  Broadband  Normal  Mode  Representation 

Normal  mode  representations  are  useful  in  describing  the  acoustic  pressure  field  in 
a  variety  of  underwater  environments.  Several  standard  textbooks  develop  acoustic 
mode  theory  in  detail  [3,  4,  5].  The  following  discussion  reviews  the  basic  concepts, 
focusing  primarily  on  modal  representations  for  broadband  signals  in  deep  ocean 
environments  such  as  those  encountered  in  ATOC.  This  section  is  split  into  two 
parts:  the  first  describes  the  use  of  the  modes  as  a  “local”  basis  set  for  the  pressure 
field;  the  second  discusses  mode  propagation  in  a  variety  of  environments. 

2.1.1  “Local”  Orthonormal  Basis 

Normal  modes  are  the  eigenfunctions  of  the  ocean  waveguide,  which  are  derived 
from  the  frequency  domain  wave  equation  (Helmholtz  equation).  At  each  frequency, 
a  mode  is  characterized  by  its  wavenumber  km  and  its  modeshape  <f)m.  For  a  given 
environment,  defined  by  the  sound  speed  profile  and  boundary  conditions,  the  modes 
satisfy  a  second-order  eigenvalue  equation,  e.g.,  in  cylindrical  coordinates  (assuming 
constant  density): 

^fa>  +  z’  n>  “  *£(r- n)]  *»(>■.  52)  =  0;  k(r,  z,  Q)  =  AL 

C  r’(2.1) 

In  Eq.  2.1,  Q  is  the  temporal  frequency,  c(r,z )  is  the  sound  speed  as  a  function 
of  range  r  and  depth  z,  and  k(r,  z,  fl)  is  the  medium  wavenumber.  The  modal 
wavenumber  (km)  determines  propagation  characteristics,  such  as  phase  and  group 
speeds,  and  the  modeshape  determines  the  spatial  distribution  of  pressure  due  to 
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each  mode.  These  shapes  are  orthogonal  functions,  scaled  such  that 

X  r^max 

~  <t>m(tt,z)<f>n(n,z)dz  =  6(m- n),  (2.2) 

p  J  0 

where  p  is  the  density  of  water. 

Since  the  modes  are  an  orthonormal  basis  for  narrowband  signals,  the  pressure 
field  at  coordinates  (r,  z)  can  be  represented  as  the  weighted  sum 

P(r,  z,  ft)  =  Y,  Mr.  Q,)<f>m{r,  z,  f2)  (2.3) 

m 

where  am  is  the  frequency-dependent  coefficient  for  mode  m.  In  general,  the  sum 
in  Eq.  2.3  is  infinite,  although  in  most  realistic  environments  only  a  finite  number 
of  modes  contribute  significantly  to  the  field.  The  remaining  “leaky”  modes  have 
complex  wavenumbers  and  suffer  exponential  losses  as  they  propagate,  thus  their 
contributions  are  negligible  in  the  far-field  of  the  source.  This  thesis  focuses  on  long- 
range  propagation  scenarios  where  it  is  reasonable  to  represent  the  pressure  field 
with  a  finite  set  of  modes. 

Time-  and  frequency-domain  representations  of  the  pressure  are  related  via  Fouri¬ 
er  synthesis,  i.e.,  the  time  series  for  a  receiver  at  range  r  and  depth  z  is 

V>(r,z,t)  =  i-  /np(r,z,n)^dn  =  ejntdU.  (2.4) 

Similarly,  the  inverse  Fourier  transform  of  the  frequency-dependent  mode  coefficient 
am(fi)  in  Eq.  2.3  defines  the  time  series  associated  with  mode  m  at  range  r: 

t)  =  ■—  J  am(r,  Sl)ejUtdVL.  (2.5) 

Limits  of  integration  in  the  above  equations  are  determined  by  the  source  bandwidth 
and  the  frequency  range  over  which  the  relevant  modes  are  propagating. 

The  following  example  illustrates  some  of  the  essential  characteristics  of  the  nor- 
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Sound  speed  profile  at  Hawaii 


Figure  2-1:  Sound  speed  at  the  ATOC  Hawaii  array  (20.2°N,  -154.0°E).  The  solid 
line  is  the  profile  obtained  using  temperature  and  salinity  from  the  Levitus  database 
for  the  winter  season.  The  dashed  line  corresponds  to  the  profile  computed  from  a 
CTD  measurement  taken  at  the  time  of  array  deployment  (mid-November  1995). 

mal  modes  using  the  environment  at  the  ATOC  Hawaii  array.  Figure  2-1  shows 
the  sound  speed  profile  for  this  location,  computed  from  Levitus  climatological  da¬ 
ta  [6,  7].  For  reference,  the  plot  also  includes  the  profile  derived  from  environmental 
measurements  taken  during  the  array  deployment.  Figure  2-2  shows  the  first  10 
modeshapes  at  75  Hz  for  both  of  these  environments.1  Note  that  each  mode  samples 
a  different  part  of  the  waveguide:  the  lowest  modes  are  concentrated  around  the 
sound  channel  axis,  while  the  higher  order  modes  cover  greater  extents  of  the  water 
column.  The  modes  of  these  two  environments  are  qualitatively  similar,  however  the 
shapes  do  reflect  the  differences  in  sound  speed  (up  to  1  m/s  near  the  channel  axis). 

Since  the  medium  wavenumber  depends  on  frequency  as  well  as  sound  speed, 
the  modes  are  functions  of  frequency,  in  general.  To  demonstrate  this,  Figure  2- 
3  compares  the  modeshapes  at  60  and  90  Hz  for  the  Levitus  environment.  The 

1  Unless  otherwise  noted,  mode  functions  are  calculated  using  Baggeroer’s  Prufer  code  [8]. 
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Speed  (m/s)  Mode  number 

Figure  2-2:  Comparison  of  Levitus  environment  and  measured  environment  at  Hawai- 
i.  The  left  panel  shows  a  closeup  of  the  sound  speed  profiles  around  the  channel  axis. 
The  right  panel  shows  the  modeshapes  for  the  first  10  modes  at  75  Hz  in  each  envi¬ 
ronment. 

frequency  range  on  the  plot  corresponds  to  the  approximate  bandwidth  of  the  ATOC 
source.  Over  this  30  Hz  interval,  mode  1  varies  slightly,  whereas  mode  10  changes 
quite  significantly.  In  general,  the  environment  and  the  source  bandwidth  determine 
the  extent  of  modal  frequency  dependence.  As  this  example  clearly  shows,  modal 
frequency  variations  are  an  important  factor  to  consider  in  the  ATOC  experiment. 

The  formulation  in  Equations  2. 1-2.5  emphasizes  that  the  modes  are  an  orthonor¬ 
mal  basis  for  a  particular  environment,  defined  by  c(r,z).  Based  on  their  spatial 
distributions,  the  lowest  modes  (eigenfunctions)  can  provide  a  compact  description 
of  acoustic  energy  concentrated  around  the  sound  channel  axis.  Beyond  being  useful 
as  a  “local”  basis  set,  the  modes  are  interesting  because  they  are  strongly  connected 
with  the  propagation  of  signals  in  the  ocean  waveguide.  The  following  section  dis¬ 
cusses  how  modes  propagate  in  a  variety  of  different  environments  and  raises  some 
questions  about  long-range  sound  transmissions. 
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Figure  2-3:  Modeshapes  for  the  first  10  modes  of  the  Hawaii-Levitus  environment  at 
60  Hz  and  90  Hz.  The  plot  shows  the  upper  2500  meters  of  the  5250  meter  waveguide. 

2.1.2  Mode  Propagation 

From  a  simple  input/output  viewpoint,  the  underwater  channel  transforms  the  mode 
signals  excited  by  the  source  into  a  modal  time  series  at  the  receiver.  Assuming  a 

finite  number  of  propagating  modes,  a  concise  frequency-domain  description  of  this 
system  is 

a[r,fi]  =  T[fi]a[0,fi]  (2.6) 

where  a[0,ft]  is  a  vector  of  mode  amplitudes  at  the  source  ( r  =  0)  and  a[r,Q] 
represents  the  corresponding  vector  at  the  receiver.  For  a  point  source,  the  modes 
are  excited  at  levels  proportional  to  the  source  spectrum  Ssic  and  the  amplitude  of 
the  modeshape  at  the  source  depth  zs  i.e., 

a[0,  tt]  =  jgyW"»(r  =  0>*»»n)  {‘l.'l) 

p{zs ) 
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In  Eq.  2.6,  T[fi]  is  a  square  matrix  that  defines  the  transformation  of  the  spectral 
amplitudes,  according  to  how  the  modes  propagate  in  the  channel.  There  are  two 
broad  classes  of  propagation  environments  to  consider:  range-independent  and  range- 
dependent. 


Range-Independent  Environments 

Given  a  fixed  sound  speed  profile  and  bottom  depth,  the  modeshapes  and  wavenum¬ 
bers  are  independent  of  range.  In  this  case,  the  modes  propagate  without  exchanging 
energy,  i.e.,  T  is  a  diagonal  matrix:2 


0 

e jkmr 

0 


(2.8) 


In  this  type  of  waveguide,  each  mode  is  a  standing  wave  in  depth  that  propagates 
outward  from  the  source  with  a  group  velocity  equal  to  d£l/dkm.  In  general,  group 
velocity  varies  with  mode  number  and  frequency,  meaning  that  the  channel  is  dis¬ 
persive.  For  a  deep  water  environment,  the  low  modes  travel  slowest,  since  they 
represent  energy  trapped  around  the  sound  speed  minimum;  higher  modes  travel 
faster.  In  a  deep  channel,  modal  group  velocity  typically  decreases  with  frequency. 


Range-Dependent  Environments 

For  realistic  ocean  waveguides,  the  environment  is  a  function  of  range,  or  more  gen¬ 
erally,  a  function  of  range  and  azimuth.  When  the  medium  is  inhomogeneous  due  to 
variations  in  sound  speed  and/or  bathymetry,  the  modes  no  longer  propagate  inde¬ 
pendently.  Instead  there  is  coupling  of  energy  among  the  modes,  meaning  that  the  T 
matrix  has  non-zero  off-diagonal  terms.  A  range-dependent  waveguide  can  be  mod- 

2Eq.  2.8  assumes  the  receiver  is  in  the  farfield  of  the  source.  See  [4,  5]  for  a  discussion  of 
range-invariant  waveguides. 
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eled  using  a  cascade  of  range-independent  segments.  In  this  type  of  model,  boundary 
conditions  at  the  segment  interfaces  determine  the  mode  coupling  coefficients. 

Since  the  coupled-mode  approach  leads  to  a  computationally-intensive  implemen¬ 
tation,  it  is  useful  to  consider  a  simplification.  The  adiabatic  approximation  assumes 
the  range  dependence  is  weak  and  neglects  the  coupling  terms,  reducing  T  to  a  di¬ 
agonal  matrix.  Under  this  assumption,  each  propagating  mode  adapts  with  range 
(changes  shape  and  wavenumber),  but  does  not  transfer  energy  into  the  other  modes. 
For  an  adiabatic  model  the  range-averaged  wavenumber, 


_  1  r 
r  Jo 


(2.9) 


determines  the  phase  and  group  speeds  for  mode  m,  thus  the  adiabatic  propagation 
matrix,  TAd,  is  simply  TRI  with  km  replaced  by  km.  Obviously  the  validity  of  the 
adiabatic  assumption  is  related  to  the  nature  of  the  inhomogeneities  in  the  medium. 
Desaubies  has  analyzed  this  problem  in  detail,  concluding  that  the  approximation’s 
accuracy  depends  strongly  on  frequency,  mode  number,  range  and  the  acoustic  quan¬ 
tity  of  interest,  e.g.,  intensity,  phase  travel  time  [9,  10]. 

In  long-range  experiments,  there  are  several  types  of  inhomogeneities  that  may 
cause  mode  coupling.  Consider  the  environment  along  the  geodesic  connecting  ATOC 
source  at  Pioneer  Seamount  to  the  receiving  array  at  Hawaii,  shown  in  Fig.  2-4.  The 
plot  shows  how  the  Levitus  (winter)  sound  speed  profiles  change  over  the  3515  km 
path.  Variability  is  concentrated  in  the  upper  500  m  of  the  water  column  and  is  a 
relatively  mild  function  of  range.  The  figure  also  displays  the  bathymetry  for  this 
section  of  the  ocean.  In  a  deep  water  environment,  changes  in  the  bottom  are  unlikely 
to  affect  the  axial  modes  since  they  do  not  interact  with  the  waveguide  boundaries. 
For  these  modes,  the  most  significant  feature  of  this  path  is  the  rapid  dropoff  in 
the  vicinity  of  Pioneer  Seamount.  The  bottom-mounted  source  (935.5  m)  does  not 
directly  excite  the  lowest  modes.  Instead  they  are  excited  by  energy  that  couples 
from  the  higher  order  modes  as  the  sound  propagates  downslope.  Chapter  4  explores 
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Avg.  SSP  Difference  from  mean  profile 


Figure  2-4:  ATOC  environment:  geodesic  path  between  source  at  Pioneer  Seamount 
(off  California)  and  the  receiving  array  near  Hawaii.  The  left  panel  is  the  average 
sound  speed  profile  over  the  path,  computed  using  235  sections  (~  15  km  apart).  The 
right  panel  shows  the  differences  between  the  mean  profile  and  the  Levitus  (winter) 
profile  for  each  of  the  sections.  Depth  of  the  ocean  bottom,  shown  in  black,  is  taken 
from  bathymetric  surveys  of  Pioneer  Seamount  [1]  and  the  ETOPO-5  topography 
database  [2]. 

the  issue  of  bathymetric  coupling  and  its  implications  for  the  mode  arrivals  measured 
in  ATOC. 

As  indicated  in  Chapter  1,  internal  waves  are  expected  to  be  the  primary  source 
of  mode  coupling  in  long-range  propagation  experiments.  In  the  deep  ocean,  in¬ 
ternal  wave  variability  is  typically  modeled  using  the  empirical  Garrett-Munk  spec¬ 
trum  [11,  12].  Figure  2-5  shows  one  realization  of  sound  speed  perturbations  due 
to  internal  wave  fluctuations,  computed  using  the  method  of  Colosi  and  Brown  [13]. 
The  calculation  assumes  the  internal  waves  are  1/2  Garrett-Munk  strength.  Note 
that  the  variability  is  greatest  in  upper  part  of  water  column. 

Before  reviewing  what  is  known  about  mode  propagation  through  random  internal 
wave  fields,  it  is  useful  to  consider  two  simulation  examples.  The  first  shows  how 
broadband  signals  propagate  through  the  slowly-range- varying  Levitus  environment. 
As  a  comparison,  the  second  example  illustrates  the  effects  of  internal  waves  by 
adding  the  sound  speed  perturbations  of  Fig.  2-5  to  the  background  environment. 


30 


Internal  wave  sound  speed  perturbations  -  single  realization 
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Figure  2-5:  Sound  speed  perturbations  due  to  internal  waves  at  1/2  Garrett-Munk 
strength 

For  simplicity,  both  simulations  use  an  axial  (rather  than  bottom-mounted)  source 
and  ignore  the  seamount. 

To  develop  some  general  intuition  about  propagation  in  the  unperturbed  Levitus 
environment,  consider  the  adiabatic  group  velocities.  Figure  Figure  2-6  shows  the  ve¬ 
locities  for  the  first  40  modes,  derived  from  the  average  wavenumbers  for  the  3515  km 
path.  On  the  plot,  the  bottom  line  on  the  plot  corresponds  to  mode  1  and  the  top 
line  corresponds  to  mode  40.  Group  velocity  decreases  with  frequency  and  increases 
as  a  function  of  mode  number,  as  is  typical  in  a  deep  water  channel.  Figure  2-7  shows 
the  predicted  spread  of  mode  arrival  times  at  3515  km  range,  assuming  a  bandlimited 
source  spectrum  between  60  and  90  Hz  and  adiabatic  propagation.  As  expected  from 
the  group  velocity  curves,  the  high  modes  arrive  first  and  exhibit  the  most  dispersion, 
while  the  low  modes  arrive  last  and  are  less  spread.  Mode  1  is  undispersed  since  its 
group  speed  is  approximately  constant  with  respect  to  frequency. 

Figure  2-8  shows  the  results  of  a  broadband  parabolic  equation  (PE)  simulation 
through  the  Levitus  background  environment.  The  top  plot  is  the  pressure  time  series 
at  the  Hawaii  array  location,  calculated  using  the  RAM  PE  code  [14]. 3  Individual 

Simulations  are  explained  in  more  detail  in  Chapter  4. 
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Range-averaged  group  velocity:  CA-HI  path 


Figure  2-6:  Range-averaged  group  velocities  for  the  first  40  modes  of  the  CA-HI 
Levitus  environment 


Adiabatic  predictions  ot  mode  spread  at  3515  km 
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Figure  2-7:  Adiabatic  predictions  of  mode  time  spread  due  to  dispersion 
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Figure  2-8:  PE  simulation  through  Levitus  environment.  The  top  panel  is  the  re¬ 
ceived  pressure  on  a  40-element  VLA;  the  bottom  panel  shows  the  corresponding 
arrival  time  series  in  the  first  10  modes 
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Figure  2-9:  PE  simulation  through  Levitus  environment  plus  internal  waves 
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mode  arrivals  are  evident  in  the  pressure  field,  e.g.,  modes  5,  3,  1.  Below  the  pressure 
plot  is  the  modal  time  series  obtained  by  projecting  the  PE  field  onto  the  functions  at 
the  receiver.  The  modes  are  obviously  arriving  in  order  from  highest  to  lowest.  The 
constructive  interference  of  the  higher  modes  result  in  the  planewave  (ray)  arrivals 
in  the  early  part  of  the  reception.4  Each  of  the  low  modes  appears  to  have  a  single 
dominant  peak,  which  arrives  at  the  predicted  adiabatic  arrival  time.  This  implies 
that  that  the  range  variations  in  the  sound  speed  do  not  result  in  mode  coupling. 
Dispersion  characteristics  of  the  waveguide  are  evident. 

In  contrast,  Fig.  2-9  shows  the  analogous  results  of  the  PE  simulation  for  the 
Levitus  environment  plus  internal  waves.  The  picture  is  quite  different.  Individual 
modes  are  no  longer  visible  in  the  pressure  time  series.  From  the  mode  time  series,  it 
appears  that  instead  of  a  single,  dispersive  arrival  in  each  mode,  there  are  multiple 
arrivals.  This  “modal  multipath”  creates  the  complicated  interference  patterns  seen 
in  the  pressure  waveforms.  Comparing  these  two  examples  to  the  ATOC  reception 
shown  in  Fig.  1-2  of  Chapter  1  reveals  that  the  experimental  measurement  more 
closely  resembles  the  internal  wave  simulation. 

From  a  theoretical  standpoint,  the  effects  of  internal  waves  on  long-range  sound 
propagation  are  not  fully  understood.  Most  previous  work  has  focused  on  the  ray 
arrivals  because  they  are  amenable  to  analysis  via  the  geometric  optics  approxima¬ 
tion.  The  monograph  by  Flatte  et  al.  summarizes  the  path  integral  theory  that 
predicts  the  fluctuations  and  coherence  of  resolved  rays  [17].  No  corresponding  the¬ 
ory  exists  for  predicting  the  behavior  of  the  mode  arrivals.  The  following  discussion 
reviews  important  results  regarding  mode  propagation  through  internal  waves  (but 
does  not  attempt  a  comprehensive  overview  of  work  in  this  area). 

In  two  seminal  papers,  Dozier  and  Tappert  derive  statistics  for  the  modal  in¬ 
tensities  of  narrowband  signals  propagating  in  a  random  ocean  [18,  19].  Based  on 
theoretical  work  and  a  set  of  numerical  simulations,  they  conclude  that  scattering 

4  A  number  of  references  explore  the  ray-mode  analogy  in  detail,  e.g.,  [15,  16,  3]. 
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eventually  results  in  an  equipartition  of  energy  among  the  modes.  To  make  the 
problem  tractable,  the  authors  rely  on  a  number  of  key  assumptions,  which  may  not 
be  valid  in  real  ocean  environments.  Notably  they  assume  that  the  acoustic  modes 
are  mutually  incoherent  (phase-random)  and  that  there  is  no  loss  of  energy  into  the 
bottom.  Regarding  the  first  assumption,  Nechaev  has  shown  that  partial  coherence 
of  the  modes  can  prevent  the  equipartition  of  energy  predicted  by  Dozier  and  Tap- 
pert  [20,  21].  Nechaev’s  analytical  results  indicate  that  decorrelation  of  neighboring 
modes  can  occur  more  slowly  than  the  randomization  of  the  overall  field  and  that 
scattered  modal  energy  can  form  a  stable  interference  structure. 

A  series  of  papers  in  the  Russian  literature  have  investigated  the  degradation  of 
mode  coherence  by  internal  waves  and  the  resulting  implications  for  various  sig¬ 
nal  processing  methods,  e.g.,  matched  filtering  [22],  horizontal  array  beamform¬ 
ing  [23,  24],  and  vertical  array  beamforming  [25].  Recently,  Sazontov  has  developed 
an  approximate  analytic  method  for  computing  the  modal  cross-coherences  and  using 
them  to  calculate  the  mutual  coherence  function  for  the  total  field  [26].  Gorodet¬ 
skaya  et  al.  provide  an  excellent  introduction  to  this  technique,  applying  it  to  a 
study  of  horizontal  and  vertical  array  gain  limitations  due  to  internal  wave  fluctu¬ 
ations  [27],  At  present,  it  is  not  known  how  well  the  approximate  expressions  for 
coherence  agree  with  experimental  data. 

Prior  to  ATOC,  there  have  been  very  few  opportunities  to  observe  the  axial 
arrivals  at  megameter  ranges.  Researchers  have  relied  heavily  on  numerical  simula¬ 
tions  to  test  theories  about  the  late-arriving  mode  energy.  In  one  of  the  first  looks  at 
experimental  data,  Colosi  et  al.  compare  pressure  measurements  from  the  1000  km 
SLICE89  experiment  to  broadband  PE  simulations  [28],  Their  results  show  that  the 
broadening  of  the  transmission  finale  in  the  data  is  attributable  to  internal  waves. 
This  smearing  in  depth  of  the  final  axial  arrivals  is  due  to  an  exchange  of  energy  a- 
mong  the  modes.  Colosi  and  Flatte  explore  the  subject  of  mode  coupling  via  internal 
waves  using  PE  simulations  designed  to  model  certain  aspects  of  the  ATOC  experi- 
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ment  [29].  They  show  that  mode  propagation  through  these  random  fields  is  strongly 
non- adiabatic  and  quantify  the  travel-time  bias/spread  and  intensity  fluctuations  for 
the  modes.  According  to  Colosi  et  aUs  recent  review  article,  internal-wave-induced 
mode  coupling,  while  definitely  an  issue  at  75  Hz,  may  be  significantly  reduced  at 
lower  frequencies,  e.g.,  28  Hz  [30]. 

Internal  waves  can  obviously  limit  the  effectiveness  of  tomographic  inversions  or 
MFP  applications  since  it  hampers  the  ability  to  associate  an  arrival  with  a  particular 
path  through  the  ocean.  As  indicated  by  this  overview,  there  is  much  left  to  learn 
about  broadband  mode  propagation  through  internal  wave  fields.  Characterizing 
the  mode  arrival  structure  is  a  prerequisite  to  using  the  modes  in  tomography  or 
source  localization.  The  ATOC  experiment  is  the  first  to  have  mode-resolving  arrays 
deployed  to  measure  axial  modes  at  megameter  ranges  over  a  period  of  months.  The 
following  section  describes  some  specific  questions  that  this  thesis  proposes  to  explore 
using  the  ATOC  receptions. 

Questions  About  Long-Range  Mode  Propagation 

This  thesis  seeks  to  address  the  following  questions  concerning  long-range  mode 
propagation.  By  answering  these  questions  we  hope  to  gain  insight  into  how  to 
identify  appropriate  observables  for  tomography  and  other  applications. 

First,  a  general  question  about  axial  arrivals  at  megameter  ranges: 

•  How  to  characterize  the  mode  arrival  structure? 

-  is  each  mode  dominated  by  a  single,  dispersive  arrival? 

-  is  there  multipath? 

-  are  the  dispersion  characteristics  of  the  channel  evident? 

•  How  do  the  mode  signals  vary  with  time? 

The  next  question  requires  a  different  approach  than  previous  researchers  have 
taken,  namely  it  requires  short-time  frequency  decompositions. 
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•  Can  the  characteristics  of  individual  rnultipaths  within  a  mode  be  measured? 
How  do  they  compare  with  the  ray  arrivals? 

-  temporal  coherence? 

-  frequency  coherence?  resolvable  multipath? 

-  fluctuation  statistics 

•  How  does  the  downslope  propagation/bottom  interaction  near  the  source  affect 
the  initial  mode  excitations.  In  turn,  how  does  that  affect  the  receptions  at 
megameter  ranges? 

For  an  experiment  like  ATOC  where  there  are  many  questions  about  the  forward 
propagation,  it  is  crucial  to  design  a  mode  estimator  that  doesn’t  assume  any  a 
priori  knowledge  of  the  arrival  structure  and  to  thoroughly  analyze  the  estimator’s 
behavior.  The  following  section  defines  the  broadband  mode  estimation  problem  and 
identifies  the  various  factors  that  determine  mode  resolution. 


2.2  Broadband  Mode  Estimation  Problem 

In  general,  the  low  order  modes  are  not  temporally  resolvable  [31],  meaning  that 
they  are  not  directly  observable  in  the  time  series  from  a  single  hydrophone.  In¬ 
stead,  vertical  arrays  can  be  used  to  separate  the  mode  signals  based  on  their  spatial 
characteristics.  This  approach  relies  heavily  on  the  orthogonality  of  the  sampled 
modes. 

Using  Eq.  2.3,  the  noisy  pressure  field  at  frequency  ft,  sampled  by  an  iV-element 
vertical  array,  may  be  written: 


p(r,  *i,  ft) 

— 

<t>i{r,  Z\ ,  ft) 

•••  <?W(r,2i,ft) 

a\  (r,  ft) 

+ 

Tl(Z\ ,  ft) 

p{r,zN,Ci ) 

(f>i  (r,  zN,  ft) 

•••  <f>M(r,zN,ty 

aM{r,  ft) 

ti(zn,  ft) 

(2.10) 
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or  in  vector  notation: 

p[r,ft]  =  $[r,f2]a[r,ft]  +  n[ft].  (2.11) 

is  the  matrix  of  sampled  modeshapes,  a  is  the  vector  of  mode  amplitudes,  and  n 
is  the  vector  of  observation  noise.  Eq.  2.10  assumes  that  the  signal  consists  of  M 
propagating  modes.  In  the  case  of  a  broadband  source,  the  array  actually  measures 
a  vector  time  series,  i.e.,  from  Eq.  2.4, 

'&(r,t)  =  f  p[r,  Q,]e3mdQ  =  f  ($[r,  fi]a[r, 0]  +  n[fi])  ejntdQ.  (2.12) 

j  vi  J  vt 

This  thesis  considers  the  problem  of  estimating  the  mode  signals  (i.e.,  a[fl])  from 
noisy  measurements  of  the  pressure  field.  A  number  of  important  signal  processing  is¬ 
sues  arise  in  designing  broadband  mode  estimators.  The  rest  of  this  section  discusses 
these  issues  in  detail,  using  the  ATOC  experiment  as  the  motivating  example. 

Modal  frequency  dependence  is  the  most  important  issue  to  consider  in  broadband 
mode  estimation.  As  Fig.  2-3  demonstrates,  modeshapes  can  vary  significantly  across 
a  30  Hz  source  band.  Clearly,  spatial  processing  must  be  done  on  a  set  of  subbands 
to  avoid  mismatch  problems.  Since  a  combination  of  temporal  and  spatial  filtering 
is  required,  there  are  time/frequency  tradeoffs  to  make.  Good  time  resolution  is 
desirable  for  resolving  the  individual  multipaths  within  a  mode.  The  allowable  widths 
of  the  subbands  is  determined  by  the  environment,  i.e.,  the  modeshape  variations. 

On  a  band-by-band  basis,  mode  estimation  reduces  to  a  classic  linear  inverse 
problem,  discussed  in  many  areas  of  the  literature,  e.g.,  estimation  theory  [32],  geo¬ 
physical  inverse  theory  [33].  In  the  narrowband  mode  estimation  problem,  the  key 
issue  to  consider  is  the  degrees  of  freedom  of  the  sampled  modeshape  matrix.  This  de¬ 
termines  how  well  the  processor  can  resolve  a  mode  and  reject  noise.  From  the  point 
of  view  of  estimating  one  mode,  there  are  two  types  of  noise  to  consider.  The  first 
component  consists  of  signals  propagating  in  the  other  modes  -  this  is  structured  in¬ 
terference  and  may  be  correlated  with  this  signal.  The  second  is  measurement  noise 
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that  uncorrelated  with  the  signal,  e.g.,  shipping  noise,  sensor  noise.  Note  that  in 
terms  of  interference  rejection,  we  assume  that  some  time  windowing  can  be  done  to 
limit  the  number  of  modes  contained  in  the  measurement,  be.,  the  ray  arrivals  (high 
order  modes)  can  be  ignored  by  time-gating. 

In  addition  to  time/frequency  tradeoffs  and  degrees  of  freedom  concerns,  two 
other  issues  arise  in  implementing  mode  estimators  for  realistic  experiments.  The 
first  concerns  arrays  that  are  not  perfectly  vertical.  This  can  be  modeled  by  using 
a  complex  modeshape  matrix  -  the  phase  terms  represent  the  timing  corrections 
required  for  each  mode.  In  planewave  beamforming,  this  is  known  as  the  array 
transit  time  problem.  It  is  important  to  quantify  the  limitations  it  places  on  the 
processor  for  modes.  Usually,  reliable  mooring  motion  estimates  are  available,  so  the 
problem  is  one  of  correcting  for  known  delays. 

A  second  practical  concern  relates  to  the  dependence  of  the  mode  on  the  local 
environment  at  the  array.  Mode  environmental  dependence  is  the  key  to  using  them 
in  tomography  or  matched  field  processing,  but  can  be  a  hindrance  when  the  receiver 
environment  is  not  exactly  known.  Uncertainty  in  our  knowledge  of  $[f2]  affects 
ability  to  resolve  the  mode  signals  accurately.  Consider  the  measured  and  archival 
profiles  shown  in  Fig.  2-2.  In  this  case,  the  archival  profile  does  not  adequately 
represent  the  modes  of  the  measured  environment.  In  ATOC  the  problem  is  that 
we  only  have  two  measurements  of  the  environment:  one  at  deployment  and  one  at 
recovery.  The  time  lapse  between  those  two  is  approximately  9  months.  It  is  expected 
that  mesoscale  effects  and  seasonal  changes  affect  the  profile  over  the  course  of  the 
experiment,  so  that  there  may  be  mismatch  between  the  modes  of  the  measured 
profile  and  the  true  profile.  It  is  important  to  quantify  the  effects  of  mismatch  on 
mode  processing. 

In  summary,  there  are  a  number  of  important  design  considerations  in  broadband 
mode  estimation.  The  following  section  reviews  relevant  literature  on  mode  filter¬ 
ing  techniques  and  their  application  to  tomography,  focusing  in  particular  on  how 
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researchers  have  approached  the  design  issues  outlined  above. 


2.3  Previous  Work 

Most  previous  research  has  focused  on  experimental  settings  where  a  narrowband 
assumption  is  valid,  meaning  that  either  the  source  is  continuous  wave  (CW)  or 
the  variations  in  the  mode  functions  across  the  source  band  are  negligible.  In  that 
context,  Ferris  and  Ingenito  describe  the  first  use  of  vertical  arrays  for  studying 
mode  propagation  in  a  real  ocean  environment  [34,  35].  Their  experiment  consists 
of  transmitting  gated  sinusoidal  pulses  in  a  shallow  water  channel  to  receivers  at 
ranges  on  the  order  of  10  km  from  the  source.  They  examine  the  mode  structure  of 
the  arriving  signals  using  both  a  temporal  filtering  method,  which  takes  advantage  of 
the  differences  in  modal  group  velocities,  and  a  spatial  filtering  method,  which  relies 
on  the  the  orthogonality  of  the  sampled  modeshapes.  For  the  latter,  the  authors  use 
a  filter  matched  to  the  modeshapes  at  the  center  frequency  of  the  transmitted  pulse. 
It  is  easy  to  show  that  the  matched  filter  (sometimes  called  the  sampled  modeshapes 
filter)  is  the  optimal  linear  mode  filter  for  detecting/estimating  any  single  mode 
in  a  background  of  spatially  white  noise  [32],  The  main  problem  with  this  filter 
is  that  it  typically  has  poor  interference  rejection  capabilities  because  it  requires 
the  orthogonality  of  the  modes  to  separate  them.  While  the  modes  are  orthogonal 
functions  over  the  continuous  aperture,  the  modeshapes,  as  sampled  by  the  array, 
may  not  be. 

A  solution  to  the  problem  of  modal  crosstalk  is  to  use  a  least  squares  approach 
to  solve  the  mode  estimation  problem.  Tindle  et  al.  are  the  first  to  introduce  this 
method  in  the  context  of  the  underwater  acoustics  [36].  In  this  case,  the  mode  filter 
is  the  pseudo-inverse  of  the  sampled  modeshape  matrix.  The  primary  advantage  of 
this  filter  is  that  the  weight  vector  for  a  particular  mode  has  nulls  in  the  directions 
of  the  other  modes  included  in  the  estimate.  Tindle  et  al.  do  note  the  connection 


40 


of  least  squares  mode  estimation  and  beamforming,  but  do  not  explore  this  idea  in 
detail. 

Polcari  explores  some  of  the  tradeoffs  of  single  beam  vs.  multiple  beam  mode 
filters  in  his  study  of  mode  coherence  in  the  Arctic  Ocean  [37].  In  addition  to 
the  least  squares  formulation,  he  considers  a  maximum  likelihood  method  that  uses 
estimates  of  the  spatial  covariance  matrix  (requiring  that  the  signal  be  stationary 
over  sufficiently  long  intervals).  Polcari  also  discusses  how  to  choose  the  number  of 
modes  to  include  in  multiple  mode  estimators. 

The  pseudo-inverse  filter  is  now  a  standard  in  the  acoustic  community,  as  indi¬ 
cated  by  its  use  in  a  number  of  matched  mode  processing  applications  [38,  39,  40]. 
The  difficulty  with  this  LS  approach  is  that  the  pseudo-inverse  may  be  poorly- 
conditioned.  This  is  a  common  problem  in  inverse  theory  [33].  A  number  of  authors 
have  considered  it  specifically  in  the  context  of  mode  processing.  Yang  suggests  using 
rank-reduced  singular  value  decompositions  to  solve  for  the  mode  amplitudes  [41]. 
Voronovich  et  al.  review  a  number  of  standard  approaches  to  solving  this  prob¬ 
lem,  but  do  not  offer  specific  recommendations  about  which  to  use  [42].  Buck  et 
al.  present  a  unified  framework  for  the  narrowband  mode  estimation  problem  [43] 
and  use  it  to  analyze  the  performance  of  several  common  mode  filters.  In  particular, 
they  discuss  the  sampled  modeshapes  filter  and  the  pseudo-inverse  filter  in  terms  of 
tradeoffs  in  interference  rejection  and  noise  rejection.  They  develop  a  maximum  a 
posteriori  mode  filter  which  gracefully  transitions  between  these  two  extremes.  Note 
that  the  MAP  approach  uses  knowledge  of  the  signal  and  noise  covariance  matri¬ 
ces;  estimating  these  quantities  typically  requires  some  assumptions  about  signal 
stationarity. 

There  are  numerous  examples  of  narrowband  mode  processing  being  applied  to 
experimental  data.  In  the  context  of  long-range  propagation,  there  are  two  ex¬ 
periments  of  note:  the  Heard  Island  Feasibility  Test  (HIFT)  and  the  Trans-Arctic 
Propagation  (TAP)  experiment.  Baggeroer  et  al.  discuss  the  use  of  least  squares 
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modal  filtering  in  the  context  of  HIFT,  which  used  a  narrowband  source  operating 
at  57  Hz  and  transmitted  across  a  17,000  km  path  to  a  vertical  array [44].  The  failure 
of  a  large  number  of  hydrophones  made  the  modal  estimation  for  HIFT  significantly 
more  difficult.  Sperry  considers  this  problem  in  detail  and  implements  a  damped 
least  squares  filter  in  order  to  perform  the  inversion  [45].  He  discusses  in  some  detail 
how  to  choose  the  number  of  modes  to  include  in  the  estimate  based  on  how  well 
the  array  samples  the  water  column.  The  recent  TAP  experiment  examined  mode 
arrivals  at  megameter  ranges  in  the  Arctic  [46].  Both  CW  and  maximal-length  se¬ 
quences  were  used.  Narrowband  processing  could  be  used  for  the  latter  because  the 
modes  did  not  vary  significantly  across  the  frequency  band  of  the  source  (centered 
at  19.6  Hz). 

In  cases  where  modal  frequency  variations  are  non-negligible,  researchers  have 
generally  approached  the  problem  by  using  an  FFT  to  separate  the  signal  into  fre¬ 
quency  bins,  doing  mode  processing  for  each  bin,  and  obtaining  a  time  series  via 
inverse  FFT.  In  one  shallow  water  experiment,  Chen  and  Lu  use  bandpass  filters 
to  facilitate  mode  processing  for  an  explosive  source  [47].  They  process  the  data  in 
bands  where  the  modeshapes  may  be  assumed  constant,  but  do  not  offer  criteria  for 
choosing  the  width  or  center  frequencies  of  these  bands.  Numerous  authors,  including 
Romm  [48],  Yang  [49],  Casey  [50]  and  Chiu  et  al.  [51],  have  implemented  broadband 
mode  estimators  using  FFT’s.  Since  they  compute  a  single  transform  for  each  receiv¬ 
er’s  time  series,  their  method  does  not  provide  a  frequency  vs.  time  decomposition 
of  the  modal  structure.  None  of  these  researchers  have  discussed  how  the  frequency 
resolution  imposed  by  the  length  of  the  FFT  affects  the  estimates.  Sutton  et  al.  have 
suggested  using  a  combination  of  time- windowing  and  Fourier  transform-based  mode 
filtering  to  resolve  modes  using  a  sparse  array  [52].  Their  article  does  not  present  the 
details  of  time  windowing  or  explore  the  implications  of  window  length  on  estimator 
resolution.  Heaney  has  used  a  time-windowing  approach  as  well  for  an  analysis  of 
the  ATOC  Engineering  test  data  [53,  54].  His  windows  for  the  FFT  are  on  the  order 
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of  3-4  seconds  long,  meaning  that  he  is  not  able  to  look  at  the  behavior  of  individual 
multipaths  in  the  Fourier  decomposition.  Heaney  uses  mode  filtering  as  a  part  of 
matched  field  inversions  for  source  location  and  internal  wave  strength.  He  does  not 
study  the  mode  arrival  structure  in  great  detail. 

Another  issue  that  needs  to  be  considered  in  broadband  mode  estimation  is  en¬ 
vironmental  mismatch.  While  mismatch  problems  in  general  have  been  thoroughly 
investigated  in  the  context  of  matched  field  processing  applications  (see  [55]  for  an 
overview),  there  has  been  little  attention  focused  on  modal  mismatch  due  to  changes 
in  the  sound  speed  profile.  Tolstoy  considers  the  effects  of  sound  speed  mismatch 
on  matched  field  localization  using  a  normal  mode  model,  but  does  not  discuss  the 
effects  on  the  individual  mode  coefficients  [56]. 

Although  Munk  and  Wunsch  first  suggested  using  modes  for  tomography  in 
1979  [57]  and  expanded  upon  the  theory  in  1983  [16],  most  experimental  applica¬ 
tions  of  tomography  have  been  based  on  ray  theory.  A  ray  is  the  result  of  the 
coherent  interference  of  a  set  of  higher  order  modes.  The  rays  form  a  stable  arrival 
structure  that  is  well-suited  to  analysis  via  planewave  beamforming  techniques.  See 
the  recent  book  by  Munk  et  al.  for  more  information  [15].  By  comparison  the  mode 
tomography  literature  is  somewhat  limited.  The  following  paragraphs  summarize 
the  relevant  publications. 

Romm  s  Ph.D.  thesis  investigates  using  modal  group  velocity  perturbations  in  to¬ 
mographic  inversions  for  sound  speed  profile  [48].  He  uses  basic  linear  inverse  theory 
to  derive  the  relationship  between  mode  group  speeds  and  sound  speed.  Simulations 
for  a  Greenland  Sea  environment  are  used  to  illustrate  the  usefulness  of  mode  to¬ 
mography.  Romm  briefly  discusses  the  required  acoustic  signal  processing,  but  his 
mode  filtering  techniques  are  quite  limited  and  poorly  explained. 

Shang  has  developed  a  full-wave  method  for  ocean  acoustic  tomography  using 
adiabatic  normal  mode  theory  and  low  frequency  sources  [58].  He  describes  how 
to  invert  for  the  sound  speed  profile  by  using  modal  phase  difference  perturbations 
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for  CW  sources  and  modal  travel  time  perturbations  for  broadband  sources.  Shang 
assumes  that  these  perturbations  can  be  measured  using  array  processing  algorithms 
or  by  using  temporal  separation  to  localize  the  pulse  arrivals  based  on  the  time  series 
from  a  single  hydrophone.  He  does  not  comment  on  how  well  the  modes  must  be 
resolved  in  order  for  the  tomographic  inversion  to  be  successful.  The  paper  includes 
a  simulation  which  demonstrates  the  feasibility  of  the  inversion  but  does  not  involve 
any  acoustic  array  processing. 

In  a  recent  set  of  papers,  Shang  and  Wang  have  applied  matched  mode  processing 
techniques  to  the  tomographic  inversion  of  an  El  Nino  profile  [59,  60].  They  base 
their  simulation  on  measured  sound  speed  profiles  from  an  actual  El  Nino  event, 
but  the  travel  time  perturbations  used  in  the  example  are  theoretical  values.  This 
research  largely  ignores  the  practical  issues  involved  in  estimating  the  mode  arrival 
times.  The  authors  do  indicate  a  way  to  compute  the  resolution  of  the  time  delay 
estimate.  Their  measure  of  resolution  comes  from  computing  the  Cramer-Rao  bound 
for  the  problem.  It  is  unclear  whether  this  bound  can  be  achieved,  with  a  practical 
system  operating  in  a  range-dependent  ocean  environment. 

The  paper  by  Sutton  et  al.  is  the  only  experiment  to  date  where  mode  time  delays 
have  been  incorporated  into  a  tomographic  inverse  [52].  They  used  a  broadband 
source  with  center  frequency  of  250  Hz  and  transmitted  across  a  106  km  path  in 
the  Greenland  Sea.  The  vertical  receiving  array  was  sparse,  consisting  of  only  6 
elements.  The  authors  use  a  least  squares  method  along  with  some  time-windowing 
to  determine  the  mode  waveforms. 

Based  on  this  review,  there  is  much  to  be  done  yet  on  broadband  mode  estimation. 
Specifically,  none  of  the  broadband  mode  processing  schemes  proposed  thus  far  have 
addressed  the  fundamental  question  of  what  frequency  resolution  is  required  for 
accurate  mode  estimation.  One  purpose  of  this  thesis  is  to  develop  and  analyze  a 
general  framework  for  broadband  mode  processing. 
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Figure  2-10:  Proposed  short-time  Fourier  mode  processing  framework 

2.4  Approach 


From  the  overview  of  previous  work,  it  is  clear  that  no  one  has  quantified  the 
time/frequency  resolution  attainable  in  mode  estimation.  To  address  these  issues, 
this  thesis  proposes  to  develop  the  short-time  Fourier  mode  processing  framework 
shown  in  Fig.  2-10.  The  basic  idea  is  to  separate  the  pressure  signals  into  a  set  of 
subbands  using  a  filterbank  and  then  do  mode  processing  in  each  subband.  The  out¬ 
put  is  a  time  series  of  mode  amplitudes  (time-varying  spectra)  that  can  be  used  to 
examine  the  frequency-dependent  structure  in  the  signals,  e.g.,  to  quantify  dispersion 
or  look  for  frequency-selective  fading. 

This  method  has  several  advantages.  First,  it  allows  us  to  examine  the  multi¬ 
ple  arrivals  within  a  single  mode  individually,  rather  than  combining  them  by  using 
a  large  FFT  for  the  required  frequency  decomposition.  This  permits  looking  for 
frequency-coherent  arrivals  within  the  multipath  arrival  pattern  and  studying  the 
fluctuation  characteristics  of  a  single  arrival  within  in  a  mode  rather  than  just  s- 
tudying  the  fluctuations  of  the  entire  mode  signal.  Note  that  the  bandpass  filtering 
effectively  limits  the  number  of  modes  that  are  present  at  each  time  step;  this  is  a 
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generalization  of  the  time-windowing  approach  proposed  by  others.  Finally,  it  is  im¬ 
portant  to  note  that  all  the  previously-described  methods  fit  within  this  framework. 

The  following  chapter  develops  the  short-time  Fourier  mode  processing  framework 
and  explores  the  temporal  and  spatial  resolution  tradeoffs.  It  also  discusses  some 
implementation  issues  related  to  non-vertical  arrays  and  environmental  uncertainty. 
Chapter  4  applies  the  short-time  Fourier  techniques  to  the  ATOC  data  set. 
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Chapter  3 


Short-time  Fourier  Mode 
Processing 


As  Chapter  2  demonstrates,  modal  variations  across  frequency  often  require  that 
spatial  processing  be  done  on  a  set  of  sub-bands  to  avoid  modeshape  mismatch 
problems.  The  purpose  of  this  chapter  is  to  present  a  framework  for  broadband 
mode  estimation  and  to  explore  the  time/frequency  resolution  tradeoffs  inherent  in 
the  processing  of  transient  mode  signals.  From  a  signal  processing  perspective,  the 
short-time  Fourier  transform  (STFT)  is  a  natural  way  to  approach  the  frequency  de¬ 
composition  required  in  mode  estimation.  The  STFT  provides  a  general  framework 
for  analyzing  the  time-  and  frequency-domain  properties  of  modal  beamforming  al¬ 
gorithms.  A  significant  advantage  to  this  approach  is  that  all  the  broadband  mode 
processing  methods  described  in  the  previous  chapter  (Section  2.3)  fit  conveniently 
within  the  basic  STFT  structure. 

The  discussion  of  short-time  Fourier  mode  processing  is  organized  as  follows. 
Section  3.1  gives  an  overview  of  short-time  Fourier  mode  processing  and  indicates 
some  of  the  important  design  tradeoffs.  Following  that,  Section  3.2  derives  equa¬ 
tions  for  the  processor  and  analyzes  performance  using  the  ATOC  experiment  as  a 
design  example.  Section  3.3  demonstrates  STFT-based  mode  processing  for  a  simple 
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adiabatic  propagation  environment.  Analysis  of  the  adiabatic  case  provides  useful 
insights  about  the  short-time  estimates.  Building  on  that  example,  Section  3.4  de¬ 
rives  the  mooring  corrections  required  for  tilted  arrays.  Section  3.5  briefly  discusses 
the  impact  of  environmental  uncertainty  on  the  STFT  mode  processor.  Finally, 
Section  3.6  summarizes  the  key  results  of  this  chapter. 

3.1  Overview  of  the  STFT  Framework 

The  STFT  is  a  standard  signal  processing  technique  for  examining  the  characteristics 
of  transient  or  time-varying  signals.  Short-time  Fourier  analysis  consists  of  comput¬ 
ing  discrete  Fourier  transforms  for  a  sequence  of  finite-length  data  segments.  There 
are  two  equally-valid  interpretations  of  the  resulting  time-dependent  spectrum:  1) 
as  the  output  of  a  filterbank,  or  2)  as  the  output  of  a  windowed  FFT  operation. 
Short-time  Fourier  synthesis  is  the  process  of  reconstructing  a  signal  from  its  time- 
varying  spectral  components.  Allen  and  Rabiner  [61]  and  Nawab  and  Quatieri  [62] 
provide  an  excellent  general  introduction  to  the  STFT;  other  authors  concentrate  on 
specific  applications.  For  example,  Rabiner  and  Schafer  [63]  discuss  the  short-time 
Fourier  method  in  the  context  of  speech  processing,  and  Johnson  and  Dudgeon  [64] 
describe  how  it  applies  to  planewave  beamforming.  In  the  latter  case,  the  STFT 
provides  a  convenient  framework  for  separating  a  multichannel  signal  into  subbands 
prior  to  spatial  filtering.  Similarly,  the  short-time  Fourier  transform  is  a  useful  way 
to  formulate  the  frequency  decomposition  required  in  broadband  mode  estimation. 

In  STFT-based  mode  analysis,  the  processor  separates  the  received  pressure  into 
a  set  of  subbands  and  computes  mode  estimates  for  each  subband.  Figure  3-1  illus¬ 
trates  these  steps  and  introduces  some  notation.  The  input  to  the  filterbank  is  \I f[l], 
a  sampled  vector  time  series  from  an  N-element  receiving  array.  As  shown,  the  filter¬ 
ing  operation  consists  of  complex  demodulation  followed  by  lowpass  filtering;  this  is 
equivalent  to  bandpass  filtering  followed  by  demodulation.  The  result  of  the  STFT 
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STFT  Analysis  Modal  Beamforming 


Figure  3-1:  Block  diagram  of  STFT-based  mode  processor 


analysis  step  is  an  7V-point  complex  vector  time  series  of  pressures  for  each  band: 
P where  k  denotes  the  band  (bin)  number  and  l  is  the  time  index.  Using  a  set 
of  narrowband  modal  beamformers,  the  processor  computes  the  time-varying  mode 
coefficients  from  the  filtered  pressures  in  each  band.  The  Af- point  vector  a.(U  [/]  is 
the  estimated  short-time  Fourier  transform  of  the  first  M  modes  in  the  bin  centered 
around  uk. 

Within  the  STFT  framework,  the  length  of  the  lowpass  filter,  HtP[u],  determines 
the  time  and  frequency  resolution  of  the  estimates.  In  the  context  of  mode  process¬ 
ing,  there  are  important  filter  length  tradeoffs  to  consider.  Long  filters  (equivalent 
to  taking  large  FFT’s)  have  good  frequency  resolution,  implying  small  operating 
bandwidths  for  the  spatial  filter  W  [uk\.  The  disadvantage  of  long  filters  is  that  they 
smear  the  arrivals  in  time.  Short  filters  provide  much  better  temporal  resolution, 
but  they  have  wider  passbands,  meaning  that  the  processor  may  be  more  sensitive 
to  the  frequency-dependent  variations  of  the  modeshapes. 
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3.2  Mode  Processing  With  Vertical  Arrays 


As  outlined  above,  the  short-time  Fourier  approach  to  mode  processing  is  conceptu¬ 
ally  simple:  separate  the  signal  into  subbands  and  estimate  the  mode  signals  in  each 
subband  using  a  narrowband  spatial  filter.  This  section  describes  the  implementa¬ 
tion  of  the  STFT  framework  for  a  perfectly  vertical  array  (i.e.,  no  tilt)  in  a  known 
environment.  A  subsequent  section  discusses  the  corrections  that  are  necessary  for 
tilted  arrays  (3.4). 

3.2.1  Broadband  Processor  Derivation 

The  purpose  of  deriving  equations  for  the  STFT  processor  is  to  facilitate  a  broadband 
performance  analysis  in  terms  of  the  time/frequency  tradeoffs  and  other  criteria 
established  in  Section  2.2.  This  derivation  uses  discrete-time  (DT)  representations, 
where  l  is  the  time  index  and  ui  is  the  DT  frequency  variable.  To  simplify  notation, 
the  mode  parameters  (shapes  and  wavenumbers)  are  written  as  functions  of  u>;  more 
precisely,  they  are  functions  of  the  corresponding  analog  frequency  12  =  u  ■  fs,  where 
fs  is  the  sample  frequency. 

Analogous  to  Eq.  2.12,  the  sampled  time  series  of  received  pressures  at  a  vertical 
array  can  be  written  using  a  discrete  Fourier  transform  (DFT)  representation,  i.e., 
the  /-point  DFT: 

1  _____  O  7T  ? 

»[!)  =  jEPNe* «i  =  -j-  >=0,...,/-l.  (3.1) 

i 

To  avoid  complications,  /  must  be  large  enough  to  represent  the  entire  reception 
(including  high  order  mode  arrivals)  and  the  subsequent  filtering  operations  without 
time  aliasing.  For  convenience,  the  following  derivation  assumes  that  the  DFT  length 
corresponds  to  the  number  of  time  samples  recorded  for  each  reception. 

The  first  step  in  processing  is  complex  demodulation,  shown  in  Fig.  3-1  as  mul¬ 
tiplication  by  e~juJkl.  After  an  appropriate  change  of  variables,  the  demodulated 
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pressure  in  the  kth.  bin  can  be  written 


^demodW  =  J  pK  +  0Ji]ejLOil.  (3.2) 

Note  that  Equation  3.2  assumes  that  the  demodulation  time  reference  is  the  first 
sample  of  the  reception  (l  =  0).  Choosing  a  different  reference  point  simply  in¬ 
troduces  an  additional  complex  exponential  term.  The  demodulation  reference  is 
important  because  the  phase  of  the  short-time  Fourier  mode  estimates  depends  on 
it. 

Following  demodulation  the  processor  lowpass-filters  each  bin,  resulting  in 

P(fc)  M  =  J  #lpMpK  +  ojt\eMl.  (3.3) 

i 

Based  on  the  development  in  Chapter  2,  the  received  pressure  at  a  particular  fre¬ 
quency  consists  of  a  finite  sum  of  modes  plus  noise,  i.e., 

PM  =  $[u/]a[w]  +  n[w].  (3.4) 

$  represents  the  matrix  of  sampled  modeshapes  at  the  receiver,  a  is  the  vector  of 
frequency-dependent  mode  amplitudes,  and  n  is  the  noise  vector. 

Substituting  this  representation  for  the  received  pressure  into  Eq.  3.3  yields 

pW[l]  =  -  J2  HLp[ui]3>[u}k  +  u;t]a[wfc  +  uji]eMl  +  j  ^  HLT[uji\n[ojk  +  Ui]ejWil.  (3.5) 

In  the  above  equation,  the  first  summation  is  the  filtered  signal  component  of  the 
pressure  field;  the  second  summation  is  the  filtered  noise  component.  Assuming  that 
the  modeshapes  are  approximately  constant  over  the  bandwidth  of  the  lowpass  filter: 

+Wj]  «  $[wfc],  (3.6) 
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then  the  measured  pressure  in  the  kth  band  becomes 


pW[/]  =  ^K]a«[/]  +  nW[/]  (3.7) 

where 

a(fc)W  =  j  HhP[uji\a[u}k  +  0Ji]eMl,  (3.8) 

n(fc)W  =  j  ^Lp[w,]n[u;jt  4-  u^e3"'1  (3.9) 

are  the  bandpass-filtered  (and  demodulated)  modal  signal  and  noise  components, 
respectively.  In  this  type  of  mode  filtering,  the  goal  is  to  estimate  the  time-varying 
mode  coefficients  a^[Z]  from  the  noisy  pressure  measurements  p^[Z]  for  each  bin. 

3.2.2  Narrowband  Mode  Filters 

Essentially,  the  STFT  approach  reduces  the  broadband  mode  estimation  problem 
to  a  set  of  narrowband  problems  where  the  measurement  takes  the  form  of  Eq.  3.7. 
Based  on  that  equation,  narrowband  mode  filtering  is  equivalent  to  determining  the 
parameters  of  a  linear  model.  This  type  of  estimation  problem  is  quite  common,  and 
there  are  a  number  of  available  techniques  for  solving  it,  e.g.,  least  squares  methods 
or  Wiener  filtering.  Selection  of  a  particular  approach  strongly  depends  on  what  is 
known  or  can  be  reasonably  assumed  about  the  signal  and  noise  characteristics.  The 
following  paragraph  states  the  assumptions  this  thesis  makes  in  designing  a  bank  of 
narrowband  mode  processors  for  the  ATOC  data. 

The  filtered  pressures  are  complex  envelopes,  resulting  from  the  demodulation 
of  real  bandpass  signals.  Both  the  mode  signals,  a^[Z],  and  the  noise,  n(fe)[/],  are 
complex  time  series.  This  thesis  assumes  that  the  noise  is  a  zero-mean  vector  random 
process,  independent  of  the  mode  signals.  Unlike  the  noise,  very  little  is  assumed 
about  the  modal  time  series.  Since  the  purpose  of  this  research  is  to  study  the  low 
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mode  arrivals  at  megameter  ranges,  it  is  important  to  avoid  biasing  those  results  by 
assuming  a  particular  temporal  or  spatial  structure  for  these  signals.  Specifically, 
this  means  that  the  processor  developed  below  does  not  rely  on  adiabatic  predictions 
of  travel  time  or  dispersion  characteristics  to  estimate  the  modes.  Furthermore, 
this  thesis  does  not  consider  time-  or  data-adaptive  processing  methods  that  require 
temporal  and/or  spatial  statistics.  Such  statistics  are  not  available  a  priori  and 
cannot  be  reliably  estimated  from  the  data  until  more  is  known  about  the  underlying 
signals,  e.g.,  their  stationarity,  etc..  Although  adaptive  algorithms  are  beyond  the 
scope  of  this  work,  they  are  reasonable  extensions  to  explore  once  more  is  known 
about  long-range  mode  propagation. 

Since  time-adaption  is  not  an  option  for  the  ATOC  data,  it  makes  sense  to  esti¬ 
mate  the  modes  independently  at  each  time  step  and  to  use  the  same  mode  filter  for 
all  time  steps.  This  thesis  focuses  on  estimates  of  the  form 

&wM  =  wfp«[j]  (3.10) 

where  W£  is  an  M  x  JV  matrix  containing  the  time-invariant  mode  filter  for  the 
kth  band.  Equation  3.10  assumes  that  the  estimated  mode  coefficients  are  linear 
combinations  of  the  received  pressures.  This  assumption  is  not  very  restrictive  since 
linear  filters  represent  a  large  subclass  of  solutions.  Most  narrowband  acoustic  mode 
filtering  algorithms  described  in  the  literature  (see  Section  2.3)  can  be  written  in 
this  form.  For  the  purpose  of  designing  the  spatial  filter,  W*,  this  thesis  assumes 
that  the  matrix  of  sampled  modeshapes  at  the  receiver  (<&)  is  known  exactly.  This  is 
equivalent  to  assuming  that  accurate  measurements  of  the  sound  speed  profile  at  the 
receiver  are  available.  Section  3.5  examines  the  impact  of  violating  this  assumption. 

In  designing  the  spatial  processing  for  the  low-modes,  this  thesis  considers  two 
standard  narrowband  mode  filters:  the  matched  (sampled  modeshapes)  filter  and 
the  pseudo-inverse  (least  squares)  filter.  There  are  a  number  of  ways  to  derive  these 
filters.  This  chapter  approaches  the  problem  from  an  array  processing  standpoint, 
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which  emphasizes  some  of  the  implicit  assumptions  being  made  about  the  underlying 
spatial  structure.  The  derivations  below  provide  a  complimentary  perspective  to  the 
estimation  theory  approach  to  mode  filtering  that  is  commonly  used  in  the  literature 
(see  [43]  for  a  summary) . 

An  important  concept  in  these  derivations  is  the  array  gain,  which  represents 
the  improvement  in  the  signal-to-noise  ratio  (SNR)  due  to  processing.  It  is  typically 
defined  as  the  ratio  of  the  SNR  at  the  output  of  a  beamformer  to  the  SNR  at  a 
single  sensor.  Since  the  signal  and  noise  characteristics  often  vary  across  an  array, 
the  input  SNR  is  taken  to  be  the  arithmetic  average  of  the  single-phone  SNR’s:1 


SNRin 


}_  A  (signal  power)n 
N  (noise  power)  n  ' 


(3.11) 


Note  that  in  the  case  of  modal  beamforming,  the  signal  levels  vary  from  one  sensor 
to  another  because  the  modeshapes  are  functions  of  depth. 

White  noise  gain  represents  the  gain  of  the  processor  when  the  noise  is  spatially 
white.  For  the  mode  processing  problem,  the  white  noise  gain  in  mode  m  is  defined 
as 


GW  =  N 


(3.12) 


where  wm  is  the  weight  vector  (filter)  for  the  mth  mode.  Application  of  the  Schwartz 
inequality  shows  that  the  maximum  value  of  the  white  noise  gain  is  N,  the  number 
of  sensors  in  the  array.  Even  if  the  noise  isn’t  spatially  white,  the  white  noise  gain 
provides  a  useful  measure  of  the  sensitivity  of  the  processor,  as  discussed  by  Cox  et 
al.  [65]. 


Matched  (Sampled  Modeshapes)  Filter 

The  matched  filter  (MF)  results  from  choosing  the  weight  vector  for  mode  m  which 
maximizes  white  noise  gain  while  maintaining  a  unit  gain  in  the  desired  mode.  Max- 

1The  geometric  average  is  used  in  some  applications. 
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imizmg  white  noise  gain  subject  to  a  unity  gain  constraint  is  mathematically  equiv¬ 
alent  to  minimizing  the  squared  length  of  the  weight  vector  subject  to  the  same  gain 
constraint,  thus  the  optimization  problem  becomes 


mm  wm  subject  to  w®0TO  =  l. 


(3.13) 


Standard  optimization  techniques  yield  the  following  solution: 


-0 


H 


m  't'mi 

<Pm<Pm 


(3.14) 


or  in  terms  of  the  weight  matrix  for  M  modes 


WH 


0 


0 


eh 


0m0m  . 


(3.15) 


where  E  is  the  sampled  modeshape  matrix  containing  the  first  M  modes,  he.,  the 
first  M  columns  of  Matched  filters  are  commonly  used  in  detection/estimation 
problems  and  their  behavior  is  well-understood.  It  is  clear  from  equations  3.12  and 
3.14  that  the  matched  filter  achieves  the  maximum  possible  white  noise  gain  for 
an  array  of  a  given  length  N .  Although  this  filter  is  optimal  in  the  sense  that  it 
maximizes  Gw,  it  does  not  explicitly  prevent  the  signal  in  one  mode  from  leaking 
into  another.  Instead,  it  must  rely  on  the  orthogonality  of  the  modes  to  separate 
them.  It  is  important  to  note  that  while  the  modeshapes  are  orthogonal  functions  of 
the  continuous  variable  z,  the  sampled  modeshape  vectors  are  not  guaranteed  to  have 
this  property.  The  degree  of  orthogonality  of  the  sampled  modeshapes  determines 
the  crosstalk  rejection  capabilities  of  the  matched  filter. 
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Pseudo-inverse  (Least  Squares)  Filter 

The  pseudo-inverse  (PI)  filter  results  from  trying  to  constrain  mode  leakage  by  plac¬ 
ing  nulls  in  the  modal  beampattern  at  the  locations  of  a  set  of  interfering  modes. 
In  this  case,  the  optimization  problem  consists  of  maximizing  the  white  noise  gain 
(minimizing  the  weight  vector  length)  subject  to  multiple  constraints,  i.e., 


minw^wm  subject  to 


Wm0m  =  1, 

-0  1  <  U  <  M. 


It  is  useful  to  rewrite  the  problem 


(3.16) 


minw^wm  subject  to  w^E  =  =  [°  0  0  0]  (3,17) 

771th  position 

where  E  contains  the  first  M  columns  of  the  sampled  modeshape  matrix  and  c  is 
an  M-point  column  vector  with  a  one  in  the  mth  position  and  zeros  everywhere 
else.  Equation  3.17  can  be  solved  using  standard  optimization  methods  involving 
LaGrange  multipliers.  Assuming  that  the  matrix  E  has  rank  M,  the  weight  vector 
for  the  mth  mode  is 

(EHE)-1EH.  (3.18) 

Equation  3.18  corresponds  to  one  row  of  the  pseudo-inverse  of  the  the  sampled  mode- 
shapes  matrix  containing  the  first  M  modes,  thus  WH  is  simply 

WH  =  (ehe)_1Eh.  (3.19) 


Provided  that  the  matrix  E,  which  contains  a  subset  of  the  sampled  modeshapes, 
has  rank  M,  the  M  —  1  null  constraints  are  met  exactly.  If  E  has  rank  less  than 
M,  a  minimum  norm  solution  of  the  constraint  equation  is  required.  The  minimum 
norm  solution  involves  the  generalized  inverse  of  E  (typically  written  in  terms  of  the 
singular  value  decomposition). 
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Another  way  to  derive  the  pseudo-inverse  mode  filter  is  to  solve  a  least  squares 
problem,  i.e.,  to  estimate  for  a  by  minimizing  the  quantity  |p-Ea|2.  This  approach  is 
the  most  common  in  the  underwater  acoustics  literature.  The  advantage  to  viewing 
mode  filtering  as  a  constrained  beamforming  problem  is  that  the  desired  spatial 
characteristics  must  be  explicitly  stated  and  that  white  noise  gain  provides  a  useful 
measure  of  the  price  paid  for  nulls  in  the  spatial  response.  White  noise  gain  of  the 
PI  filter  for  the  mth  mode  is 


N  1 

'  Cm(EHE)“1C; 


(3.20) 


From  this  expression,  it  is  clear  that  the  white  noise  gain  of  the  PI  filter  is  intimately 
connected  with  the  conditioning  of  the  pseudo-inverse.  The  term  second  term  is 
bounded  by  the  squared  singular  values  of  the  E  matrix: 


^min  ^  cT(EHE)_1C  ~  °max  (3.21) 

where  crmin  and  crmax  are  the  minimum  and  maximum  singular  values  of  E,  respec¬ 
tively.  The  white  noise  gain  of  the  PI  filter  less  than  or  equal  to  that  of  the  matched 
filter  (which  achieves  the  maximum  possible  Gw). 

Conditioning  problems  are  well  known  in  both  least  squares  problems  and  in  adap¬ 
tive  filtering.  Menke  discusses  methods  of  dealing  with  poorly-conditioned  pseudo- 
inverses  in  the  context  of  inverse  theory  [33].  Cox  provides  a  nice  discussion  of 
robustness  issues  in  the  context  of  adaptive  beamforming  [65].  There  are  many  ways 
of  dealing  with  this  problem.  Cox  provides  some  time-adaptive  approaches  (which 
are  not  appropriate  for  the  ATOC  data  due  to  the  transient  nature  of  the  signals). 
Non-adaptive  solutions  to  conditioning  problems  typically  involving  deleting  small 
singular  values  from  the  inverse,  or  added  a  small  weight  to  the  diagonal  terms  of  EHE 
before  inverting  it  in  order  to  stabilize  the  inversion.  Both  these  methods  bias  the 
results  and  should  be  used  with  care  (especially  when  not  much  is  known/understood 
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about  the  arrival  structure). 

Thus  far,  this  section  has  derived  equations  for  the  STFT  processor,  showing  how 
it  reduces  the  broadband  mode  estimation  problem  to  a  set  of  narrowband  problems. 
Two  approaches  to  narrowband  filtering  have  been  discussed.  The  following  subsec¬ 
tions  explore  narrowband  and  broadband  performance  issues,  using  a  design  example 
based  on  the  ATOC  experiment. 

ATOC  Design  Example 

Numerical  results  presented  below  use  the  Levitus  environment  for  the  California- 
Hawaii  path,  described  in  Chapter  2.  Note  that  although  the  archival  sound  speed 
differs  from  the  measured  profile  at  the  array  (see  Fig.  2-1),  these  differences  do 
not  affect  the  basic  conclusions  regarding  the  number  of  modes  that  can  be  reliably 
estimated,  the  allowable  bandwidth,  etc.  The  array  in  the  example  corresponds  to  the 
ATOC  VLA  configuration:  40  elements,  with  a  spacing  of  35  m,  spanning  the  water 
column  between  330  m  and  1695  m.  For  reference,  Fig.  3-2  shows  the  hydrophone 
locations  relative  to  the  first  20  modes  at  75  Hz  for  the  Hawaii-Levitus  environment. 
75  Hz  is  the  center  frequency  for  the  ATOC  source,  which  has  approximately  a  30  Hz 
bandwidth  (from  60-90  Hz).  The  sample  rate  for  the  experiment  is  300  Hz. 

3.2.3  Narrowband  Performance  Analysis 

The  main  purpose  of  this  section  is  to  compare  the  matched  filter  and  the  pseudo¬ 
inverse  narrowband  beamformers,  focusing  primarily  on  the  tradeoffs  in  noise  and 
interference  rejection.  Beampatterns  provide  a  useful  measure  of  the  crosstalk  rejec¬ 
tion  capabilities  of  a  spatial  processor.  In  the  context  of  modal  beamforming,  the 
beampattern  is  defined  as  201og10(WH$),  where  W  is  the  multidimensional  mode 
filter  and  <I>  is  the  matrix  of  sampled  modeshapes.  The  mth  row  of  the  beampattern 
matrix  corresponds  to  the  projection  of  the  modes  into  the  estimate  for  mode  m. 
For  the  matched  filter,  the  beampattern  corresponds  to  a  normalized  version  of 
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Figure  3-2:  Modeshapes  (75  Hz)  and  receiver  locations  (+’s)  for  the  design  example 

the  sampled  modeshape  correlation  matrix,  Figure  3-3  shows  MF  beampattern 
for  the  first  40  modes,  using  the  ATOC  array  at  Hawaii.  This  beamformer  has 
excellent  crosstalk  rejection  for  the  lowest  modes,  up  until  approximately  mode  8. 
Above  that,  the  modeshapes  sampled  by  the  array  are  obviously  correlated.  As 
a  result,  energy  from  one  mode  will  leak  into  its  neighbors.  Performance  of  this 
beamformer  degrades  significantly  for  modes  above  10;  this  is  not  surprising  since 
the  array  is  designed  to  sample  the  first  10  modes. 

Figure  3-4  shows  the  beampatterns  for  three  different  pseudo-inverse  filters,  de¬ 
signed  for  10  modes,  15  modes,  and  20  modes,  respectively.  Note  that  the  Fig.  3-4(c) 
has  a  different  dB  scale  than  figures  3-4(a)  and  3-4(b).  By  design,  the  PI  filter 
guarantees  nulls  for  a  specified  set  of  modes,  thus  the  beampatterns  have  a  diagonal 
structure.  It  is  important  to  remember  that  while  the  PI  filter  guarantees  a  fixed  set 
of  nulls,  it  does  nothing  to  prevent  higher  order  modes  not  included  in  the  pseudo¬ 
inverse  from  leaking  into  the  lower  order  modes.  The  beampatterns  show  how  the 
higher  order  modes  (up  to  40)  not  included  in  the  pseudo-inverse  alias  into  the  lower 
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Input  Mode 


Figure  3-3:  MF  Beampattern  (75  Hz)  for  the  40-element  ATOC  VLA  in  the  Hawaii- 
Levitus  environment 


modes. 

The  price  paid  for  the  nulls  in  the  PI  beampattern  is  a  loss  in  white  noise  gain. 
Figure  3-5  shows  the  white  noise  gain  as  a  function  of  mode  number  for  the  10-, 
15-,  and  20-mode  filters.  The  dotted  line  shows  the  maximum  white  noise  gain, 
101og10(Af),  which  is  achieved  by  the  matched  filter.  For  the  10-mode  PI  filter,  the 
white  noise  gain  is  equivalent  to  the  matched  filter  results.  The  15-  and  20-mode 
PI  filters  show  the  loss  in  white  noise  gain  associated  with  poor  conditioning  of  the 
pseudo-inverse.  For  the  20-mode  filter,  the  white  noise  gain  is  negative  for  some 
modes  meaning  that  the  noise  will  be  significantly  amplified  compared  to  the  signal. 

The  amount  of  crosstalk  is  also  governed  by  the  conditioning  of  the  pseudoinverse, 
i.e.,  crosstalk  increases  as  the  conditioning  degrades.  This  is  important  because, 
while  it  may  be  possible  to  assume  that  the  time-windowing  imposed  by  the  lowpass 
filter  operation  in  the  STFT  processor  eliminates  the  highest  modes  (earliest  ray 
arrivals),  the  potential  for  internal- wave-induced  mode  scattering  means  that  there 
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Input  Mode 

(a)  10-mode  PI  filter 


Input  Mode 

(b)  15-mode  PI  filter 


(c)  20-mode  PI  filter;  note  the  different  dB  scale. 


Figure  3-4:  Comparison  of  pseudo-inverse  filter  beampatterns  (75  Hz)  for  the  ATOC 
array  in  the  Hawaii-Levitus  environment 
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Mode  number 

Figure  3-5:  White  noise  gain  at  75  Hz  for  the  ATOC  array  in  the  Hawaii-Levitus 
environment 

may  be  a  large  number  of  modes  within  the  time  window  spanned  by  the  temporal 
filter.  As  Fig  3-4(c)  shows,  the  crosstalk  from  higher  order  modes  dominates  the 
estimate  for  the  20-mode  filter,  as  a  result  of  poor  matrix  conditioning.  Figure  3-6 
is  a  plot  of  the  maximum  crosstalk  as  a  function  of  the  number  of  modes  considered. 
The  x-value  represents  the  highest  mode  considered  and  the  y-value  represents  the 
peak  crosstalk  into  that  set  of  modes.  For  example,  the  leftmost  point  on  the  plot 
represents  how  much  crosstalk  there  is  into  mode  1  from  all  the  other  modes  up 
to  mode  40.  This  plot  illustrates  how  using  the  PI  strategy  of  placing  nulls  at 
neighboring  modes  in  order  to  eliminate  crosstalk  is  a  losing  game  when  there  are  a 
large  number  of  modes  propagating  in  the  waveguide.  As  the  number  of  constraints 
increases,  so  does  the  length  of  the  weight  vector,  meaning  that  modes  still  not 
included  in  the  estimate  are  amplified  even  more. 

It  is  interesting  to  note  that  the  degradation  in  conditioning  of  the  pseudo-inverse 
exhibited  in  the  20-mode  filter  is  not  significant  enough  for  a  matrix  inversion  pro¬ 
gram  (such  as  the  one  used  by  Matlab)  to  give  a  warning  about  conditioning.  From 
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Figure  3-6:  Maximum  crosstalk 


a  numerical  perspective,  the  inverse  can  still  be  computed  without  errors,  but  from 
an  array  processing  perspective,  performance  is  significantly  degraded.  If  there  were 
no  noise  in  the  reception,  then  the  mode  estimates  could  be  computed  effectively, 
even  with  this  type  of  poorly  conditioned  pseudo-inverse.  Obviously,  assuming  the 
absence  of  noise  is  unrealistic  for  any  experiment.  There  are  ways  of  mitigating 
conditioning  problems,  but  only  at  the  expense  of  biasing  the  estimates. 

This  section  has  shown  how  modeshape  orthogonality  affects  the  matched  filter 
results,  how  the  PI  filter  null  placement  strategy  can  eliminate  some  of  the  crosstalk 
problems  but  at  the  expense  of  robustness.  It  has  been  shown  that  noise  and  crosstalk 
from  higher  order  modes  both  increase  as  the  pseudo-inverse  becomes  ill-conditioned. 
Based  on  the  numerical  examples,  it  is  clear  that  the  PI  filter  for  10  modes  gives  noise 
gain  results  equal  to  that  of  the  matched  filter  while  providing  the  added  benefit  of 
eliminating  the  crosstalk  among  modes  8  through  10  which  is  present  in  the  matched 
filter.  In  terms  of  the  ATOC  experiment,  the  PI-10  filter  is  a  reasonable  choice 
for  the  best  narrowband  performance.  The  following  section  explores  broadband 
performance  issues. 
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3.2.4  Broadband  Performance  Analysis 

In  designing  a  short-time  Fourier  mode  processor,  it  is  crucial  to  have  a  measure 
of  how  well  a  narrowband  mode  filter  designed  for  frequency  u>  performs  on  the 
modes  at  a  neighboring  frequency,  u  ±  Auj.  Essentially,  the  operating  bandwidth 
determines  the  required  frequency  resolution  of  the  lowpass  filter,  which  in  turn 
defines  the  temporal  resolution  of  the  STFT  processor.  This  section  explores  these 
broadband  performance  issues.  First,  it  examines  the  frequency  response  of  both  the 
MF  and  PI  mode  filters,  using  the  ATOC  design  example.  Based  on  these  results,  it 
is  possible  to  select  an  appropriate  lowpass  filter/data  window  for  the  STFT.  Using 
that  selection,  we  discuss  the  time  resolution  of  the  STFT  processor  for  the  ATOC 
example.  Finally,  this  section  considers  the  noise  response  as  a  function  of  frequency. 
This  is  important  because  it  gives  a  measure  of  the  conditioning  of  the  mode  filter 
across  the  source  band  (allows  us  to  explore  how  the  conditioning  of  the  mode  filter 
is  affected  by  the  changes  in  the  modeshapes  across  frequency). 

Frequency  Resolution 

Consider  the  signal  component  of  the  mode  estimate  (combining  equations  3.10 
and  3.5): 

a  {k)[l]  =  W  k(jY,  #lpN<I>K  +  ^]ak  +  Wi]^  (3.22) 

Based  on  this  equation,  the  frequency  response  of  the  STFT  processor  is  primarily 
governed  by  the  matrix  W^4»[c«;fc  +  u;,],  which  represents  the  frequency-dependent 
beampattern  of  the  mode  filter  W*.  Figures  3-7  and  3-8  show  plots  of  these  beam- 
patterns  for  the  matched  filter  and  the  pseudo-inverse  filter.  Both  filters  are  designed 
using  the  Hawaii-Levitus  environment  modeshapes  at  75  Hz;  the  PI.  filter  is  designed 
for  10  modes.  Results  are  shown  for  the  30  Hz  (±15  Hz)  band  around  the  center 
frequency.  In  these  figures,  each  subplot  corresponds  to  the  beampattern  for  a  single 
mode.  The  solid  line  represents  the  response  in  the  desired  mode  and  the  dashed 
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Mode  2 


60  70  80  90  60  70  80  90 

Mode  7  Mode  8 


60  70  80  90  60  70  80  90 

Mode  9  Mode  10 


Frequency  (Hz)  Frequency  (Hz) 

Figure  3-8:  Frequency  response  of  the  PI  filter  designed  for  10  modes.  Solid  lines 
indicate  the  response  in  the  desired  mode;  dashed  lines  indicate  crosstalk  from  neigh¬ 
boring  modes  (up  to  10). 
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lines  represent  the  crosstalk  from  neighboring  modes  (from  1  to  10)  into  the  desired 
mode. 

First  consider  the  matched  filter  results  for  mode  1.  The  beamformer  has  a  flat 
response  for  this  mode  across  the  30  Hz  band.  There  are  nulls  in  the  direction  of  the 
other  modes  at  the  center  frequency,  but  crosstalk  increases  sharply  as  w  deviates 
from  the  design  frequency  (75  Hz  in  this  case).  Plots  for  the  higher  modes  indicate 
that  frequency  mismatch  can  significantly  affect  the  gain  in  the  desired  mode  as 
well  as  the  crosstalk  rejection,  e.g.,  the  gain  in  mode  10  is  down  by  10  dB  at  60 
Hz.  As  expected  from  the  narrowband  analysis,  Fig.  3-7  also  shows  that  the  MF 
beamformer  does  not  prevent  crosstalk  at  the  center  frequency  in  modes  8  through 
10.  The  beampatterns  for  the  PI  mode  filter  in  Fig.  3-8  show  similar  behavior.  Each 
filter  is  constrained  at  75  Hz  to  have  unity  gain  in  the  desired  mode  and  nulls  at  the 
other  modes,  but  these  constraints  are  not  guaranteed  to  be  maintained  for  other 
frequencies. 

Both  the  MF  and  PI  beampatterns  indicate  that  crosstalk  rejection  is  more  signif¬ 
icantly  affected  by  frequency  mismatch  than  is  the  gain  in  the  desired  mode.  Further 
insight  into  this  observation  is  obtained  by  considering  a  perturbation  theory  analysis 
of  the  modeshapes  as  a  function  of  frequency.  Appendix  A  discusses  the  application 
of  perturbation  theory  to  the  mode  eigenvalue  problem  in  general  and  presents  the 
frequency  perturbation  analysis  upon  which  the  results  below  depend.  The  basic 
idea  is  that  the  modeshape  at  frequency  w  +  Aw  can  be  represented  as  the  sum  of 
the  modeshape  at  frequency  w  plus  first-  and  second-order  correction  terms: 

</>m[w  +  Aw]  «  0m[w]  +  (~^j  (3.23) 

where  $  is  the  matrix  of  sampled  modeshapes  at  frequency  w  and  gP>  and  g$  are 
column  vectors  containing  the  first-  and  second-order  perturbation  coefficients.  As 
Eq.  3.23  indicates,  the  correction  terms  are  written  in  terms  of  the  modes  themselves, 
since  they  form  a  complete  orthonormal  basis.  In  general,  the  perturbation  of  mode 
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m  due  to  a  change  in  frequency  is  well-represented  by  a  small  set  of  mode  rrC s  nearest 
neighbors,  thus  only  a  few  coefficients  in  the  g  vectors  are  important. 

Using  the  perturbation  expansion,  first  consider  the  response  in  the  desired  mode 
i.e., 

Wm4UW  +  A^]  ~  wm<t>m M  +  (~)  w£$g^}  +  w“$gg)  (3.24) 

where  the  dependence  of  the  weight  vector  on  co  has  been  suppressed.  Assuming  that 
the  weight  vector  effectively  passes  only  the  desired  mode  at  the  center  frequency  (a 
reasonable  assumption  for  a  well-conditioned  PI  filter),  then  this  expression  reduces 
to 

+  Au]  w  1  +  g£l  (3.25) 

where  g$m  and  g$m  represent  the  amount  of  <j>m  contained  in  the  first-  and  second- 
order  perturbations  of  <f>m.  As  derived  in  Appendix  A,  g$m  is  equal  to  0,  and  g$m 
is  negative.  Thus  the  loss  in  the  gain  of  the  desired  mode  is  a  second-order  effect. 
Similarly,  looking  at  the  crosstalk  factor  shows  that 

+  Aw]  «  0  +  gfl  (3.26) 

In  this  case,  g^  is  not  generally  equal  to  zero,  thus  the  crosstalk  terms  go  as 

Perturbation  theory  results  can  be  used  to  obtain  approximate  bandwidths  for 
mode  filters.  In  practice,  it  is  often  just  as  easy  to  plot  the  response  and  choose  the 
bandwidth  based  on  that  rather  than  an  approximation.  Figure  3-9  is  a  closeup  of 
the  frequency  response  of  the  10-mode  PI  filter.  From  this  figure,  it  is  clear  that  the 
crosstalk  stays  20  dB  down  for  frequencies  within  2.5  Hz  of  the  center  frequency.  This 
seems  like  a  reasonable  3  dB  bandwidth  to  choose  for  the  lowpass  filter.  Based  on 
the  perturbation  theory  analysis,  the  crosstalk  sidelobes  increase  at  a  rate  of  20  dB 
per  decade,  thus  it  is  useful  to  select  a  lowpass  filter  whose  sidelobes  fall  off  faster 
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Frequency  (Hz) 

Figure  3-9:  Frequency  response  of  the  PI  filter  for  mode  10  (designed  using  10  modes) 
than  20  dB  per  decade. 

Although  #LP  can  be  any  filter  with  a  lowpass  characteristic,  short-time  Fourier 
analysis  typically  employs  the  same  set  of  windows  that  is  used  in  spectral  esti¬ 
mation  [66].  Harris  provides  an  extensive  list  of  window  functions  in  his  classic 
paper  [67],  As  argued  above,  the  ATOC  data  requires  a  filter  with  approximately  a 
5  Hz  (±2.5  Hz)  bandwidth  with  sidelobes  that  fall  off  at  a  rate  greater  than  20  d- 
B  per  decade.  Consider  the  Hanning  window,  which  has  sidelobes  that  fall  off  at 
60  dB  per  decade  (in  terms  of  the  power  metric,  201og10  |tfLP|).  Figure  3-10  shows 
the  impulse  and  frequency  responses  for  Hanning  windows  of  three  different  lengths, 
60  points,  120  points,  and  240  points,  assuming  the  ATOC  sample  rate  of  300  Hz. 
From  this  plot  it  appears  that  an  0.4  second  (120-pt)  Hanning  window  gives  the 
desired  bandwidth  and  sidelobe  characteristics.  It  provides  adequate  suppression  of 
mode  crosstalk  for  modes  up  to  10.  This  choice  of  filter  defines  the  time  resolution 
of  the  processor.  The  following  section  discusses  the  temporal  response  of  the  STFT 
processor. 
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Hanning  windows 


(a)  Impulse  responses 


Hanning  window  frequency  responses 


(b)  Frequency  response  magnitudes 


Figure  3-10:  Time  and  frequency  responses  for  three  Hanning  windows,  assuming  a 
300  Hz  sample  rate.  Amplitude  differences  among  the  impulse  responses  are  a  result 
of  constraining  the  response  to  be  equal  to  1  at  0  Hz. 


Time  Resolution 

To  understand  the  time  resolution  of  the  bandpass-filtered  mode  estimates,  it  is  use¬ 
ful  to  determine  the  minimum  separation  between  two  arrivals  such  that  they  are 
resolvable  as  individual  peaks  at  the  output  of  the  filter.  For  simplicity,  suppose  the 
input  signal  consists  of  two  equal-amplitude  impulses,  then  the  output  is  the  sum 
of  two  filter  impulse  responses.  Figure  3-11  depicts  the  result  of  this  superposition 
for  two  different  assumptions  about  the  spacing  of  the  impulses.  A  reasonable  cri¬ 
terion  for  resolving  the  arrivals  is  that  the  response  at  the  midpoint  between  them 
is  smaller  than  the  response  at  the  peak  locations.  This  is  true  provided  that  the 
input  separation  is  greater  than  the  “3  dB  timewidth”  of  the  window,  i.e.,  the  width 
between  the  points  where  the  window  amplitude  is  equal  to  1  /2  of  its  peak.  For  a 
Hanning  window,  the  3  dB  timewidth  is  equal  to  L/(2fs)  where  L  is  the  number 
of  points  and  fs  is  the  sampling  frequency.  This  definition  of  temporal  resolution  is 
obviously  oversimplified,  in  particular  because  the  effects  of  noise  have  been  ignored. 
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(a)  Unresolvable  (b)  Resolvable 

Figure  3-11.  Simple  illustration  of  time  resolution  imposed  by  the  lowpass  filter 

Nevertheless,  it  is  useful  as  an  approximate  measure  of  the  temporal  smearing  intro¬ 
duced  by  the  lowpass  filter.  The  taper  of  the  0.4  second  (120  point)  Hanning  window 
is  such  that  it  is  not  possible  to  resolve  arrivals  spaced  more  closely  than  0.2  seconds. 

In  addition  to  smearing  the  signals  in  time,  the  lowpass-filtering  operation  also 
introduces  a  delay.  The  group  delay  of  the  filterbank  specifies  the  time  adjustment 
required  to  accurately  relate  the  bandpass-filtered  mode  estimates  to  the  received 
pressure  time  series.  Provided  that  the  filters  have  approximately  constant  delay 
across  the  passband,  the  correction  involves  a  simple  shift  of  the  time  axis.  For¬ 
tunately,  many  of  the  standard  analysis  windows  (including  the  Hanning  window) 
qualify  as  generalized  linear  phase  systems2 .  Such  systems  are  known  to  have  a  con¬ 
stant  group  delay,  which  is  equal  to  (^f±)  fs  (where  L  is  the  filter  length).  All  of 
the  STFT  estimates  discussed  in  this  thesis  implicitly  incorporate  the  filter  delay 
correction. 

Noise  Response 

A  previous  section  has  shown  that  the  pseudo-inverse  filter  designed  for  10  modes 
is  well-conditioned  and  achieves  nearly  the  optimal  white  noise  gain.  Since  the 
modeshapes  are  functions  of  frequency,  it  is  important  to  examine  how  the  noise 


Specifically,  these  windows  are  real  and  symmetric:  h[l]  =  h[L-l-l\,  0  <  l  <  L— 1.  Oppenheim 
and  Schafer  discuss  the  properties  of  generalized  linear  phase  systems  in  detail  [68]. 


Mode  1 


Mode  10 


(a)  Mode  1 


(b)  Mode  10 


Figure  3-12:  White  noise  gain  as  a  function  of  frequency  for  the  matched  filter  and 
three  different  pseudo-inverse  filters  (for  10  modes,  12  modes,  and  15  modes). 


characteristics  vary  over  the  source  band.  Figure  3-12  shows  the  white  noise  gain  for 
modes  1  and  10,  plotted  over  the  30  Hz  interval  corresponding  to  the  ATOC  source. 
The  plots  have  4  curves:  one  for  the  matched  filter  and  three  for  different  realizations 
of  the  pseudo-inverse  filter  which  are  designed  for  10  modes,  12  modes,  and  15  modes, 
respectively.  Mode  1  is  clearly  unaffected  by  frequency,  but  conditioning  problems 
are  clearly  evident  for  mode  10  PI  filters.  It  is  useful  to  compare  these  results  with  the 
plot  of  the  modeshapes  at  60  and  90  Hz,  shown  in  Fig.  2-3  of  Chapter  2.  The  changes 
in  mode  1  with  respect  to  frequency  are  negligible  and  that  results  in  constant  white 
noise  gain  across  the  band.  Mode  10,  on  the  other  hand,  changes  significantly.  At  the 
lower  frequencies  the  modeshape  covers  a  greater  extent  of  the  water  column  than 
it  does  at  the  higher  frequencies.  As  a  result  the  array  does  not  sample  the  mode  as 
well,  therefore  the  pseudo-inverse  is  more  poorly  conditioned.  Based  on  these  results, 
the  10-mode  PI  filter  still  appears  to  be  reasonable  choice  for  the  ATOC  data  since 
its  noise  response  does  not  change  significantly  across  the  source  band. 
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Conclusions  of  the  ATOC  Design  Example 

Based  on  the  numerical  results  presented  above,  the  pseudo-inverse  filter  for  10  modes 
provides  the  best  possible  combination  of  noise  gain  and  interference  rejection.  This 
filter  has  a  bandwidth  of  approximately  5  Hz  (±2.5  Hz).  An  0.4  second  Hanning 
window  has  the  desired  bandwidth,  and  has  the  additional  advantage  of  sidelobes 
that  fall  off  faster  than  the  crosstalk  increases.  These  are  the  default  parameters  for 
the  STFT  processor  designed  for  the  ATOC  data.  Using  this  design,  the  next  section 
applies  STFT  techniques  to  analyze  a  simple  propagation  example. 

3.3  Adiabatic  Example 

The  purpose  of  this  section  is  to  illustrate  important  characteristics  of  the  STFT 
processor  using  an  adiabatic  propagation  example.  The  discussion  of  this  example 
consists  of  two  parts.  First,  Section  3.3.1  shows  the  received  pressure  and  corre¬ 
sponding  mode  estimates  for  adiabatic  propagation  through  the  California-Hawaii 
Levitus  environment.  Following  that,  Section  3.3.2  derives  analytic  expressions  for 
these  estimates  to  verify  the  numerical  results  and  to  provide  more  intuition  about 
this  type  of  mode  processing. 

3.3.1  Processing  Results 

Table  3.1  summarizes  the  parameters  for  the  adiabatic  simulation.  The  environment 
contains  a  broadband  point  source  at  700  m  depth,  located  3515.2  km  away  from  a 
40-element  receiving  array  (35  meter  spacing).  In  this  example,  the  source  transmits 
a  single  windowed-sinusoidal  pulse  with  approximately  30  Hz  bandwidth,  at  a  center 
frequency  of  75  Hz.  The  time  series  for  the  receivers  is  synthesized  using  the  first 
40  modes.  Figure  3-13  shows  the  received  pressure  as  a  function  of  time  and  depth 
for  the  adiabatic  simulation.  For  this  case,  the  modes  are  dispersed  enough  (at 
3515.2  km  range)  that  it  is  possible  to  identify  individual  modes  in  the  pressure  time 
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Figure  3-13:  STFT  example:  processing  of  adiabatic  propagation  data.  Top  plot  is 
the  received  pressure  on  a  40-element  array  as  a  function  of  time  and  depth.  The 
bottom  plots  show  the  frequency-stacked  outputs  for  modes  1  and  5,  respectively. 


Source 

depth: 

center  frequency: 
pulse: 

pulse  duration: 

700  m 

75  Hz 

triangular-windowed  sinusoid 
.11  secs  (~  30  Hz  bandwidth) 

Receiver 

range: 

3515.2  km 

no.  of  receivers: 

40 

element  spacing: 

35  m 

span: 

330  m  -  1695  m 

sample  rate: 

300  Hz 

Modes 

number  synthesized  at  receiver: 

40 

Table  3.1:  Simulation  parameters  for  the  California-Hawaii  adiabatic  propagation 
example 
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series,  e.g.,  mode  1  is  the  strong  final  arrival,  and  mode  5  is  associated  with  the  five 
strong  peaks  lined  up  right  after  2374.5  seconds. 

Mode  estimates  were  calculated  from  this  pressure  time  series  using  a  10-mode 
PI  filter  and  the  0.4  second  Hanning  window  filter  described  in  the  previous  section. 
Note  that  “critical  sampling”  for  this  Hanning  window  requires  computing  bins  every 
2.5  Hz.  Critical  sampling  refers  to  the  bin-spacing  required  in  theory  to  be  able  to 
resynthesize  a  time  series  from  the  STFT  decomposition.  For  this  example  the  STFT 
processor  used  a  bin-spacing  of  1.25  Hz;  the  additional  samples  improve  the  look  of 
the  plots,  but  do  not  improve  the  frequency  resolution.  The  bottom  two  plots  in 
Fig.  3-13  show  the  frequency-stacked  outputs  for  modes  1  and  5.  There  are  several 
interesting  features  to  note.  First,  in  this  example  these  modes  are  actually  evident  in 
the  pressure  time  series.  The  dispersion  characteristics  of  the  two  modes  are  evident 
in  the  plots.  Mode  1  is  effectively  undispersed,  i.e.,  all  it  arrives  at  the  same  time 
in  all  frequency  bins.  Mode  5  is  clearly  somewhat  dispersed  since  the  signal  at  lower 
frequencies  arrives  before  the  signal  at  higher  frequencies. 

Figure  3-14  shows  the  results  of  STFT  processing  for  the  first  10  modes.  The 
significantly  lower  signal  levels  for  modes  4,  6,  and  9  are  the  result  of  low  excitation 
of  the  mode  by  the  source  due  to  the  low  amplitude  of  these  modeshapes  at  the 
source  depth.  Modes  9  and  10  contain  some  energy  arriving  before  their  adiabatic 
arrival  times;  this  is  due  to  crosstalk  from  higher  order  modes  not  nulled  out  by  the 
10-mode  PI  filter.  As  noted  above,  the  dispersion  characteristics  of  the  modes  are 
evident  in  these  plots.  Note  that  because  of  the  temporal  smearing  caused  by  the 
filter,  the  dispersion  has  to  be  on  the  order  of  0.2  seconds  before  it  will  be  evident 
in  the  STFT  output. 

3.3.2  Analysis 

The  starting  point  for  the  analysis  of  the  adiabatic  example  is  the  filtered  pressure 
time  series  (Eq.  3.5).  Assuming  that  the  signal  consists  of  only  the  mth  mode  and 


75 


Frequency  (Hz)  Frequency  (Hz)  Frequency  (Hz)  Frequency  (Hz)  Frequency  (Hz) 


Mode  1 


Mode  2 


Mode  3  Mode  4 


Figure  3-14:  Estimated  frequency-stacked  outputs  for  the  first  10  modes,  computed 
using  a  pseudo-inverse  filter 
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that  there  is  no  additive  noise,  the  received  pressure  on  the  nth  hydrophone  in  the 
kth  band  is 


Pn]il\  -  J  J2  #LP [Wi\<l>m[Zn,Uk  +  Ui]am[u}k  +  (3.27) 

i 

As  discussed  previously,  the  lowpass  filter  has  a  generalized  linear  phase,  meaning 
that  its  frequency  response  is  of  the  form: 


#lpM  =  H[u]eju(‘L21)  =  (3.28) 

where  H[u]  is  a  real  and  even  function  of  w,  L  is  the  number  of  points  in  the 
filter,  and  tu  represents  the  delay  (in  seconds)  associated  with  filter,  i.e.  f#  =  -1-. 
For  simplicity,  this  derivation  assumes  that  the  lowpass  filter  is  designed  to  prevent 
significant  spatial  mismatch,  thus  <j)m[uk  +  Wj]  ss  for  all  tUj  passed  by  the  filter. 
Using  that  assumption  and  substituting  in  the  generalized  linear  phase  representation 
of  the  lowpass  filter  yields 


Pn\l]  ~  <t>m[zn,uk]  (^j  ^ H[u>i]am[uk  +  Wj]e  /.) eMl  (3.29) 


The  adiabatic  model  presented  in  Chapter  2  gives  an  expression  for  the  frequency- 

dependent  mode  amplitude  am.  Based  on  equations  2.6  and  2.7,  this  amplitude 
is 


^771  ^ 
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(3.30) 


In  Eq.  3.30,  the  first  term  represents  the  excitation  of  mode  m  by  a  point  source 
located  at  the  range  origin  and  the  depth  zt.  The  second  term  represents  the  effects 
of  the  channel  on  the  signal  propagating  in  mode  m.  In  this  term,  km  is  the  adiabatic 
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(range-averaged)  wavenumber,  which  is  a  function  of  frequency.3  The  final  term 
represents  the  phase  shift  required  to  time-advance  the  signal  so  that  the  7-term 
Fourier  expansion  corresponds  to  a  reception  starting  at  time  i0  seconds  instead  of 
0  seconds.  Before  substituting  the  expression  for  am  into  Eq.  3.29,  it  is  convenient 
to  define  the  following  spectral  amplitude  for  mode  m: 


Am[u),r 


Sarc[(«#m[r  =  0,Za,W 
p[zs]yJSTTkm[uj\r 


(3.31) 


For  simplicity,  assume  that  the  source  is  zerophase  and  that  its  magnitude  varies 
slowly  with  respect  to  u.4  Since  the  modeshape  and  wavenumber  are  also  smooth 
functions  of  u,  the  amplitude  Am  is  slowly-varying  function  of  frequency.  Provided 
that  the  lowpass  filter  is  sufficiently  narrowband,  the  term  Am[uk  +  w*],  required  by 
Eq.  3.29,  is  approximately  equal  to  Am[<jjk\.  With  these  assumptions,  the  pressure 
time  series  may  be  written 


P(n][l}  =  Am[uk,rn}ei”/4<t>m[zn,uk}  0)  ^ 

(3.32) 

where  rn  is  the  range  to  the  nth  hydrophone. 

The  next  step  in  the  analysis  requires  a  model  for  the  variations  of  the  modal 
wavenumber  with  frequency.  Since  the  range-averaged  wavenumber  is  not  an  analytic 
function  in  general,  it  is  necessary  to  use  a  Taylor  series  for  km.  Expanding  the 


3The  use  of  an  overbar  (km)  to  indicate  the  range-averaged  wavenumber  has  been  dropped  to 
simplify  notation. 

4  This  is  equivalent  to  assuming  either  that  the  source  transmits  a  symmetric  pulse  or  that 
the  receiver  uses  pulse  compression  (a  temporal  matched  filter  applied  to  the  time  series  on  each 
hydrophone).  In  the  latter  case  Ssrc  represents  the  Fourier  transform  of  an  autocorrelation  function, 
which  is  symmetric.  To  be  more  precise,  the  source  has  a  generalized  linear  phase  representation 
(to  account  for  time  shifts  of  the  pulse).  There  is  no  loss  of  generality  in  assuming  zerophase  since 
any  linear  phase  term  can  easily  be  incorporated  into  the  term  in  Eq.  3.30. 
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Figure  3-15:  Difference  in  range-averaged  wavenumber  and  3-term  Taylor  series  ap¬ 
proximation  for  the  CA-HI  Levitus  environment.  Note  that  the  phase  differences  are 
scaled  by  the  range  (3515  km)  of  the  ATOC  experiment  and  normalized  by  ir. 


wavenumber  around  the  frequency  Dc  (continuous-time)  yields 


km[^]  ~  km[f2c]  +(D  —  Qc) 

Qc  Spin 


dkm 

dQ, 

n=nc 

+ 


(3.33) 


n=nc 


The  terms  spm,  sgm,  and  dm  represent  the  phase  slowness,  group  slowness,  and  dis¬ 
persion  coefficient,  respectively,  for  mode  m.  Qc  is  defined  to  be  the  center  frequency 
of  the  source,  e.g.,  fc  =  75  Hz  for  the  ATOC  experiment.  Figure  3-15  illustrates 
the  accuracy  of  this  approximation  for  the  CA-HI  Levitus  environment.  The  plot 
shows  the  differences,  scaled  by  range,  in  the  actual  wavenumber  for  the  path  and 
the  3-term  Taylor  series  approximation.  As  the  figure  shows,  phase  errors  due  to 
the  approximation  are  less  than  7t/4  for  the  20  Hz  band  (65-85  Hz)  for  the  first  10 
modes.  Differences  increase  substantially  outside  of  that  band,  especially  for  the 
higher  order  modes. 
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Using  the  Taylor  series  expansion,  the  second-order  approximation  to  the  kmr 
term  in  Eq.  3.32  is  (after  converting  to  DT  frequency) 

km[uk  +  Wy ]r  «  LUcspmrfs  +  (Au:k  +  uJi)sgmrfs  +  ^(Aw*  +  coi)2dmrf]  (3.34) 

where 

A  uk  —  Dk  UJC.  (3.35) 

Note  that  spmr  is  the  phase  delay  (seconds)  and  sgmr  is  the  group  delay  (seconds) 
associated  with  propagation  through  the  channel.  Substituting  Eq.  3.34  into  Eq.  3.32 
and  taking  the  terms  not  involving  Ui  outside  the  summation 

pW[l)  =  Ame^<j>m[zn, uk]  0)  £ (3.36) 

where  Bm  is  a  phase  shift: 

7T  1 

Bm[u)k,  Tn]  —  —  ^c^pmTn  )fs  A u)k(sgmTn  to)/s  —(A ujk)  dmrnfs ,  (3.37) 

tm  corresponds  to  the  arrival  time  (in  seconds)  of  the  mth  mode  at  frequency  cjk: 

Ui]  =  SgmTn  d"  A U)kdmTnf s,  (3.38) 

and  /?  is  a  dispersion  factor: 

PVn)  =  dmrnfl  (3.39) 

Equation  3.36  looks  complicated,  but  has  a  rather  simple  interpretation.  AmejBm  is  a 
complex  scale  factor  and  <f>m  is  the  spatial  weighting  associated  with  the  modeshape. 
The  summation  represents  the  inverse  Fourier  transform  of  a  pulse: 

AputaW  =  (j)e  £M  (3.40) 

*  pulse  shape  pulse  location  distortion 


80 


(3.41) 


~  h  [l  —  (tm  +  tH  -  to)fs]  =  h[l  -  lm\ 

The  first  factor  in  the  sum  determines  the  pulse  shape,  the  second  term  determines  its 
location,  and  the  third  produces  a  distortion  of  the  pulse  shape  (due  to  the  dispersion 
of  mode  m).  Since  the  passband  of  the  filter  is  assumed  to  be  very  narrow,5  the 
distortion  of  the  pulse  shape  is  minimal.  Neglecting  the  dispersion  factor 
the  pulse  is  simply  a  shifted  version  of  the  impulse  response  of  the  filter.  In  other 
words,  the  pulse  looks  like  the  window  associated  with  the  short-time  transform,  as 
is  indicated  in  the  approximate  expression  in  Eq.  3.41.  It  is  important  to  note  that  it 
is  only  the  effect  of  dispersion  on  the  pulse  shape  that  is  being  neglected  here.  Modal 
dispersion  also  determines  the  phase  of  a  mode  arrival  across  the  frequency  bins  of 
the  STFT,  and  in  general,  the  dispersion  term  in  Eq.  3.37  cannot  be  neglected. 
Substituting  Eq.  3.41  into  Eq.  3.36  yields 

Pn][l]  =  <t>m[zn}Am{rn}eJBm[rn]h  [l  -  lm[rn) } .  (3.42) 

where  the  dependence  of  (j)m,  Am,  Bm,  and  lm  on  u>k  is  implicit  and  the  dependence  on 
the  receiver  coordinates  (r„,  zn)  has  been  emphasized.  Recall  that  spatial  processing 
consists  of  computing  a  weighted  sum  of  the  pressure  over  N  sensors,  i.e., 

M  =  E  Wm[Zn]p{nk)[l}  (3.43) 

n= 1 

where  wm{zn }  is  the  nth  component  of  the  weight  vector  for  mode  m.  Substituting 
Eq.  3.42  into  Eq.  3.43  results  in  an  approximate  expression  for  the  estimated  time 
series  for  mode  m.  The  resulting  equation  illustrates  what  happens  for  both  vertical 
and  non-vertical  array  geometries. 

First,  consider  the  case  of  a  vertical  array,  where  rn  =  r  for  all  n.  Most  of  the 

5Note  that  if  the  narrowband  assumption  for  H[ui)  is  violated,  spatial  crosstalk,  rather  than 
phase  distortion,  is  likely  to  be  the  dominant  source  of  errors. 
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Figure  3-16:  Magnitude  and  phase  across  frequency  at  the  peak  arrival  time  in  the 
75  Hz  bin  for  modes  1  and  5  of  the  adiabatic  example 


terms  can  be  pulled  outside  the  sum,  therefore 

Omty]  =  Ame?Bmh[l  -  lm]  Wm[zn}<l>m\Zn}  (3.44) 

n= 1 

V - - - ' 

=  1 

The  summation  is  equal  to  one  because  of  the  unity  gain  constraint  imposed  on  the 
weight  vector.  From  Eq.  3.44,  the  mode  m  signal  in  the  kth  band  is  a  pulse  centered 
around  sample  lm.  Converting  from  samples  to  seconds,  mode  m  arrives  at  time 
t  —  tm  +  tf{  (assuming  the  source  emits  a  pulse  at  t  =  0).  The  amplitude  and  phase 
(Am  and  Bm,  respectively)  do  not  vary  over  the  extent  of  the  arrival,  i.e.,  they  are 
not  functions  of  l.  Amplitude,  phase  and  arrival  time  do  depend  on  the  frequency 
ojk  as  indicated  by  equations  3.31,  3.37,  and  3.38. 

Consider  taking  a  slice  of  the  STFT  estimate  across  frequency  at  the  peak  arrival 
time.  This  gives  a  measure  of  the  amplitude  Am  and  the  phase  Bm.  Figure  3-16 
shows  the  frequency  slices  taken  at  the  peak  arrival  times  for  modes  1  and  5  of  the 
adiabatic  example.  The  amplitude  reflects  the  spectrum  of  the  source  (the  triangular- 


82 


windowed  sinusoid).  For  mode  5,  the  amplitude  is  also  affected  by  the  fact  that  the 
peak  arrivals  in  the  bins  do  not  line  up  exactly  (due  to  dispersion  within  mode  5). 
The  latter  effect  is  minimal. 

Note  that  the  phase  of  the  arrival  across  frequency  bins  has  a  linear  term  associ¬ 
ated  with  its  arrival  time  relative  to  the  start  of  the  reception.  This  is  expected  since 
a  delay  in  the  time  domain  corresponds  to  a  phase  shift  of  the  Fourier  transform. 
For  the  results  shown  in  Fig.  3-16(b)  this  trend  has  been  removed.  In  theory,  the 
linear  component  can  be  completely  removed,  given  knowledge  of  the  arrival  time 
relative  to  the  start  of  the  processing  interval.  In  practice,  a  small  linear  component 
may  remain  due  to  the  fact  that  the  picked  peak  may  not  correspond  to  the  exact 
arrival  time  of  the  pulse  (can  be  off  due  to  noise,  etc.). 

As  long  as  dispersive  spreading  within  mode  m  is  not  greater  than  the  window 
width,  the  amount  of  dispersion  can  be  measured  using  a  frequency  slice  taken  at 
a  constant  time.  If  the  dispersive  spreading  is  larger  than  the  window  width  (0.4 
seconds  in  this  case),  then  it  would  be  necessary  to  look  along  a  slanted  line  in 
time-frequency  space.  For  the  adiabatic  modes  of  the  California-Hawaii  Levitus 
environment,  the  maximum  dispersive  spreading  up  to  mode  10  is  approximately  0.4 
seconds. 

How  do  the  adiabatic  results  generalize?  It  is  worth  pointing  out  that  even  in 
non-adiabatic  cases,  we  can  still  associate  a  magnitude,  phase,  and  time  with  an 
arrival.  In  the  general  case,  the  effective  wavenumber  would  be  a  weighted  average 
of  the  wavenumbers  of  the  modes  that  the  signal  has  traversed  the  path  in.  Removal 
of  the  linear  term  in  the  phase  is  still  possible  because  that  simply  requires  knowledge 
of  the  arrival  time  of  the  pulse.  In  this  way,  it  may  be  possible  to  get  insight  about 
the  dispersion  characteristics  of  a  particular  multipath  arrival. 
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3.4  Mooring  Corrections 


The  previous  discussion  has  assumed  a  perfectly  vertical  array  for  the  pressure  mea¬ 
surements,  i.e ,  all  sensors  located  at  the  same  range  from  the  source.  This  section 
relaxes  that  assumption  and  considers  the  practical  problem  of  compensating  for 
array  tilt.  Since  STFT  processing  involves  a  joint  time/frequency  representation  of 
the  signals,  the  mooring  corrections  involve  both  a  phase  and  a  time  adjustment. 
For  the  purposes  of  this  chapter,  it  is  assumed  that  accurate  location  information  is 
available  for  all  sensors.  This  is  a  reasonable  assumption  in  many  experiments  like 
ATOC. 

Building  upon  the  derivation  in  the  previous  section,  assume  that  the  array  is  not 
vertical.  In  cylindrical  coordinates,  this  means  that  the  distance  to  the  nth  sensor 
is  rn  =  r  -  A rn,  where  r  is  a  reference  distance.  In  this  case  the  mode  estimate 
becomes  (from  Eq.  3.43): 

a&’M  =  £  Mr  -  a -  A/m))  (3.45) 

71  —  1 

where  A Bm  and  A lm  depend  on  A rn: 

Ai?m  ^c^pmfs  T  AwfcSsmArn/s  -f-  —  (Atajt)  dmArnfs ,  (3.46) 

Alm  —  ( Sgm  T  A U7fcrfm)  fs Arn.  (3.47) 

In  reducing  this  equation,  first  consider  the  amplitude  Am,  given  by  Eq.  3.31.  Pro¬ 
vided  that  r  is  large  compared  to  A rn,  the  amplitude  is  unaffected  by  the  fact  that 
the  array  is  not  vertical.  This  is  certainly  true  in  ATOC  where  r  is  on  the  order 
of  megameters  and  A rn  on  the  order  of  100  meters.  The  phase  shift  and  time  shift 
terms,  A Bm  and  A lm  are  much  more  significant.  Note  that  the  phase  shift  A Bm  is 
simply  the  second  order  approximation  to  kmArn.  Using  this  information,  Eq.  3.45 
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becomes 


=  4eiBm  E  hll  -  Vm  -  Mm)}wl[zn]4>m[zn}e^^.  (3.48) 

n=l 

From  this  equation,  it  is  clear  that  the  array  tilt  results  in  a  phase  shift  and  a  time 
shift  shift  of  the  pulse  associated  with  each  sensor.  In  general,  compensation  for  both 
of  these  effects  is  required  for  proper  interpretation  of  the  STFT  mode  estimates. 

Based  on  Eq.  3.48,  the  way  to  implement  the  phase  adjustment  is  to  design  the 
mode  filters  for  complex  modeshapes,  i.e., 


(f>m(zn)ejkmri 

<t>m{zn)ejkmT2 

</>m(2n)ejfcmrjv 


(3.49) 


Using  complex  modeshapes  is  the  standard  way  to  incorporate  the  phase  compen¬ 
sation  into  the  MF,  PI,  or  other  mode  filter  designs.  Note  that  the  derivations 
presented  in  Section  3.2.1  are  equally- valid  for  the  case  of  complex  modeshapes. 

The  corrections  for  the  time  shift  due  to  array  tilt  are  easily  implemented  by 
adding  the  appropriate  delay  to  the  lowpass-filtering  operation.  It  is  useful  to  con¬ 
sider  whether  the  time-corrections  need  to  be  done  on  a  mode-by-mode  basis  or  not. 
The  dominant  term  in  the  time  shift,  lm  is  sgmArn  (the  dispersion  being  a  2nd  order 
effect).  Figure  3-17  shows  the  difference  in  the  delay  term  associated  with  mode  m 
and  the  delay  term  associated  with  mode  1,  i.e.,  tm  -  h  as  a  function  of  the  change 
in  range  A r„.  The  group  slownesses  used  for  this  calculation  are  those  associated 
with  the  range-averaged  Levitus  wavenumbers.  The  differences  are  less  than  3  mil¬ 
liseconds  up  to  differential  ranges  of  10  km.  For  most  vertical  array  geometries,  the 
maximum  value  of  Arn,  will  be  much  less  than  the  maximum  value  shown  here.  The 
maximum  A rn  for  the  ATOC  experiment  is  less  than  300  m.  Clearly,  for  that  con- 
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Figure  3-17:  Difference  in  delay 


figuration,  the  difference  due  to  mode  number  is  negligible.  Making  the  assumption 
that  the  time  shift  is  mode-independent  simplifies  the  processing,  since  it  means  that 
the  filterbank  will  be  identical  for  all  modes. 


3.5  Impact  of  Environmental  Uncertainty 

The  purpose  of  this  section  is  to  briefly  address  the  issue  of  environmental  uncertainty 
in  the  context  of  broadband  mode  processing,  specifically  in  terms  of  the  ATOC 
experiment.  The  modeshapes  used  in  designing  the  mode  filters  depend  on  the  local 
environment  at  the  receiving  array.  For  the  ATOC  experiment,  there  were  only  two 
measurements  of  that  environment:  one  when  the  array  was  deployed  in  November 
of  1995  and  one  when  the  array  was  recovered  in  August  of  1996.  Figure  3-18  shows 
the  sound  speed  profiles  associated  with  these  two  environmental  measurements. 
The  differences  in  these  two  measured  profiles  are  slight,  compared  to  the  difference 
between  the  measured  profiles  and  their  archival  counterparts  (see  the  comparison  of 
Levitus  and  measured  in  Chapter  2.  Still,  these  differences  do  affect  the  mode  shapes 
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Figure  3-18:  Comparison  of  deployment  and  recovery  sound  speed  profiles 

and  wavenumbers,  which  means  they  affect  the  performance  of  the  spatial  filter  in 
two  ways:  in  the  spatial  beampattern  and  in  the  mooring  corrections.  Figure  3-19 
shows  the  modeshapes  at  75  Hz  for  the  two  different  profiles.  First  consider  how 
mismatch  affects  the  beampattern.  Figure  3-20  shows  the  effective  PI  beampattern 
when  the  mode  filters  are  designed  using  the  deployment-profile  modes,  but  the  true 
underlying  modes  of  the  signal  correspond  to  those  for  the  recovery  profile.  Figure  3- 
21  shows  similar  results  for  the  matched  filter  case.  As  both  of  these  beampatterns 
illustrate,  mismatch  produces  nearest-neighbor  coupling  of  energy.  In  the  worst  case, 
neighboring  modes  are  only  about  10  dB  down  (in  power)  from  the  peak  in  the  desired 
mode.  While  there  are  slight  differences  in  the  mismatched  beampatterns  for  the  PI 
and  MF  designs,  the  magnitude  of  the  coupling  is  the  same  for  both. 

In  addition  to  altering  the  beampattern,  mismatch  may  also  affect  the  mooring 
corrections.  For  ATOC,  the  maximum  difference  in  the  wavenumbers  at  75  Hz  be¬ 
tween  the  deployment  and  recovery  profiles  is  approximately  2*  10~5,  thus  the  phase 
errors  involved  in  correcting  for  a  maximum  A rn  of  300  m  are  negligible. 

Since  only  2  measurements  of  the  receiver  environment  are  available  for  the  ATOC 
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Modeshapes  @75.00  Hz 


Figure  3-19:  Comparison  of  deployment  and  recovery  modeshapes  at  75  Hz 


experiment,  it  is  hard  to  characterize  the  maximum  possible  mismatch  (the  pro¬ 
files  during  when  the  receptions  occurred  could  be  very  different  from  either  of  the 
measurements).  The  archival  profiles  are  offset  from  both  the  measured  profiles 
by  enough  that  they  are  not  useful  for  estimating  the  profile  in  between  the  two 
measurements.  Accurate  characterization  of  the  local  environment  is  obviously  an 
important  issue  in  mode  estimation.  It  is  possible  to  envision  schemes  for  mitigating 
the  impact  of  this  mismatch,  but  they  require  more  information  and  are  beyond  the 
scope  of  this  thesis. 


3.6  Summary 

This  chapter  has  presented  a  short-time  Fourier  framework  for  broadband  mode 
processing.  The  narrowband  mode  filters  that  are  a  part  of  the  STFT  structure  were 
derived  using  a  different  approach  than  the  standard  one  in  the  acoustics  literature. 
This  approach  made  the  connection  of  mode  processing  to  optimal  beamforming. 
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Pseudo-inverse  Filter  Beampattern:  10  modes 
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Figure  3-20:  Beampattern  for  the  mismatch  case:  PI  filter  designed  with  deployment 
profile  modeshapes;  recovery  profile  modes  are  the  input 


Input  Mode 


Figure  3-21:  Beampattern  for  the  mismatch  case:  matched  filter  designed  with  de¬ 
ployment  profile  modeshapes;  recovery  profile  modes  are  the  input 
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A  thorough  performance  analysis  of  the  system  has  explored  both  narrowband  and 
broadband  issues.  Based  on  this  analysis,  the  STFT  processor  parameters  for  the 
ATOC  data  were  chosen. 
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Chapter  4 


STFT  Mode  Analysis  of  ATOC 
Receptions 


As  stated  in  Chapter  1,  the  motivation  for  developing  a  framework  for  broadband 
mode  estimation  is  to  study  the  mode  arrivals  at  megameter  ranges.  This  chapter 
presents  an  analysis  of  the  ATOC  receptions  on  the  Hawaii  VLA,  using  the  short-time 
Fourier  techniques  described  in  the  previous  chapter.  Processing  of  the  ATOC  data 
reveals  a  complicated  multipath  arrival  pattern  for  the  first  10  modes  at  3515  km 
range.  This  chapter  investigates  the  frequency  characteristics  of  the  measured  ar¬ 
rival  structure,  compares  experimental  estimates  to  simulation  results,  and  identifies 
several  useful  statistics  for  characterizing  the  low-mode  signals. 

The  outline  of  this  chapter  is  as  follows.  Section  4.1  begins  by  introducing  the 
ATOC  data  set,  focusing  on  the  receptions  at  Hawaii  considered  in  this  thesis.  Fol¬ 
lowing  that,  Section  4.2  describes  a  set  of  PE  simulations  used  for  comparison  pur¬ 
poses.  Section  4.3  provides  the  first  look  at  the  time- varying  mode  spectra  for  the 
ATOC  receptions  and  compares  the  experimental  results  to  the  STFT  estimates  for 
simulated  receptions.  Next,  Section  4.4  addresses  the  issue  of  temporal  variability  by 
comparing  ATOC  data  with  simulated  receptions,  and  computing  average  coherence 
as  a  function  of  time.  Section  4.5  examines  some  statistics  of  the  mode  estimates, 
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obtained  by  averaging  over  multiple  transmissions.  Finally,  Section  4.6  summarizes 
major  conclusions  of  this  chapter. 

4.1  ATOC  Data  Set 

As  described  in  Chapter  1,  the  ATOC  experiment  involved  a  bottom-mounted  source 
at  Pioneer  Seamount  (off  California)  broadcasting  to  two  vertical  arrays  and  a  set 
of  bottom-mounted  receivers  in  the  Northeastern  Pacific.  The  source  transmitted 
phase-encoded  pseudo-random  sequences  at  a  center  frequency  of  75  Hz,  with  a 
3  dB  bandwidth  of  37.5  Hz.  Each  transmission  consisted  of  44  repetitions  of  the 
27.28  second  pseudo-random  sequence,  corresponding  to  a  transmission  length  of 
approximately  20  minutes.  The  source  and  the  two  VLA’s  were  deployed  in  the  fall 
of  1995,  and  regular  transmissions  began  in  December  of  that  year.  Over  the  course 
of  the  experiment,  the  source  transmitted  signals  every  four  hours  during  periods  set 
by  the  Marine  Mammal  Research  Program  associated  with  ATOC. 

This  thesis  considers  only  the  data  recorded  on  the  Hawaii  VLA,  which  was 
located  at  a  range  of  3515.2  km  from  the  source.  Figure  4-1  shows  the  schedule  of 
reception  times  for  this  array  from  December  28,  1995  (yearday  362)  to  May  23, 
1996  (yearday  509). 1  Although  the  array  was  not  recovered  until  August  1996  the 
deepest  20  hydrophones  failed  sometime  after  day  509.  Since  mode  processing  is 
much  more  difficult  with  only  the  shallow  half  of  the  array  (due  to  poor  sampling  of 
the  modeshapes),  this  study  is  limited  to  188  receptions  recorded  on  the  full  array. 

For  each  transmission,  the  40-element  array  recorded  10  four-period  averages 
(over  18.2  minutes)  of  the  pseudo-random  sequence  received  on  each  hydrophone.2 
Subsequently,  the  time  series  for  each  sensor  was  demodulated  and  matched-filtered 

^his  plot  only  shows  the  times  of  188  “good”  receptions.  There  were  229  transmissions  between 
yeardays  362-509,  but  41  of  these  have  been  ignored  due  to  incomplete  or  corrupted  data. 

2  Although  the  source  transmitted  44  periods  of  the  M-sequence,  the  array  only  recorded  40  pe¬ 
riods.  The  start-time  for  recording  was  chosen  so  that  sampling  began  approximately  two  periods 
after  the  start  of  the  reception  [69]. 
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Hawaii  Reception  Schedule 
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Figure  4-1:  ATOC  transmission  schedule  through  yearday  509.  Crosses  mark  the 
time  of  each  good  reception;  receptions  with  bad  channels  have  been  eliminated  from 
the  data  set.  The  line  of  numbers  below  the  crosses  indicates  how  the  receptions  are 
divided  into  13  subgroups  for  post-processing. 

to  achieve  pulse  compression.  For  the  purposes  of  all  the  numerical  results  presented 
in  this  thesis,  the  received  pressure  time  series  consists  of  these  matched-filtered 
demodulates. 

During  the  experiment,  the  position  of  the  VLA  was  tracked  using  a  long-baseline 
acoustic  navigation  system  consisting  of  four  transponders  deployed  on  the  bottom 
and  several  interrogator  hydrophones  on  the  array.  Navigation  data  was  recorded 
immediately  before  and  after  each  reception.  This  thesis  relies  on  the  mooring  motion 
calculations  done  by  the  group  at  Scripps.  Array  positions  for  each  of  the  10  four- 
period  averages  were  determined  by  interpolating  between  the  beginning  and  ending 
locations  of  the  array. 

Before  analyzing  the  mode  arrivals  in  the  ATOC  receptions,  it  is  useful  to  have 
estimates  of  the  input  noise  levels.  The  intent  is  to  provide  some  approximate  noise 
statistics  for  use  in  setting  thresholds  for  plotting  and  peak  detection,  rather  than 
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to  undertake  a  thorough  analysis  of  the  noise  at  the  ATOC  arrays.  Since  no  sepa¬ 
rate  noise  measurements  were  made  for  the  ATOC  experiment,  noise  levels  must  be 
estimated  from  the  receptions  themselves.  The  approach  taken  here  is  to  use  two 
“noise-only”  segments  of  the  27.28  second  received  time  series  to  calculate  statistic- 
s.  For  the  Hawaii  array,  at  least  7  seconds  at  the  beginning  of  the  reception  and 
4  seconds  at  the  end  of  the  reception  do  not  appear  to  contain  any  strong  signal 
components.  These  two  sections  of  data  are  used  in  the  noise  analysis. 

Figure  4-2  shows  estimated  noise  spectra  on  two  hydrophones  for  the  first  recep¬ 
tion  on  yearday  363  (one  of  the  4-period  averages).  These  estimates  were  computed 

Hydrophone  10  Hydrophone  30 


Figure  4-2:  Estimated  noise  spectra  for  hydrophones  10  and  30.  The  solid  line 
is  computed  using  the  noise-only  section  prior  to  the  signal  arrival;  dashed  line  is 
computed  using  the  noise-only  section  after  the  signal  arrivals. 

using  the  modified  periodogram  method  of  Welch  [70,  66].  The  data  window  used  in 
the  processing  is  an  0.4  second  Hanning  window,  thus  the  frequency  resolution  of  the 
noise  estimates  is  comparable  to  that  of  the  STFT  mode  processor  (approximately 
±2.5  Hz  bandwidth).  Window  overlap  is  50%,  resulting  in  30  and  24  segments  for 
the  periodogram  averages  of  the  first  and  second  noise-only  sections,  respectively.  As 
the  figure  shows,  noise  power  estimates  for  the  first  and  second  sections  (beginning 
and  end  of  the  received  time  series,  respectively)  are  consistent.  Noise  is  obviously 
a  function  of  frequency,  with  the  levels  at  90  Hz  being  down  from  those  at  60  Hz  by 
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Figure  4-3:  Estimated  spatial  covariance  in  the  75  Hz  bin;  calculated  from  300 
snapshots.  Units  are  dB,  referenced  to  peak  value. 

approximately  6  dB. 

Based  on  Fig.  4-2,  the  noise  levels  on  hydrophone  numbers  10  and  30  are  com¬ 
parable.  To  get  a  better  idea  about  the  spatial  characteristics  of  the  noise,  consider 
the  spatial  covariance  matrix  shown  in  Fig.  4-3.  This  covariance  matrix  is  calculated 
for  the  75  Hz  bin  of  the  same  reception.  The  beginning  (30  segments)  were  used 
and  the  results  were  averaged  over  the  10  periods  (4-period  averages)  recorded  by 
the  receiver.  The  covariance  is  dominated  by  the  diagonal  terms,  with  some  nearest 
neighbor  correlation,  but  no  other  obvious  structure. 

The  results  shown  in  Figures  4-2  and  4-3  are  representative  of  the  results  for 
other  receptions  and  support  the  conclusion  that  the  noise  field  is  approximately 
spatially  white.  Since  the  results  for  the  first  noise-only  section  are  similar  to  the 
second  section,  it  is  reasonable  to  average  those  results  together.  In  addition,  the 
noise  estimates  can  be  averaged  over  hydrophones  (assuming  spatially  white)  and 
over  the  10  4-period  averages  in  order  to  obtain  a  single  average  noise  spectrum 
for  each  reception.  No  frequency-averaging  has  been  done.  Figure  4-4  shows  how 


95 


these  averaged  noise  statistics  vary  over  the  course  of  the  experiment.  For  reference, 


Peak  signal  vs.  noise  power 


Figure  4-4:  Input  noise  levels  at  60,  75,  and  90  Hz  as  a  function  of  reception  number. 

the  peak  signal  power  in  each  reception  (he.,  the  maximum  pressure  received  on  a 
single  hydrophone)  is  also  shown.  Note  that  the  noise  varies  significantly  enough 
from  reception  to  reception  that  it  is  important  to  set  plotting/detection  thresholds 
individually  for  each. 

Assuming  a  spatially  white  model  for  the  input  noise,  it  is  easy  to  compute  the 
noise  levels  at  the  output  of  each  bin  of  the  STFT  mode  processor.  If  the  gain  of 
the  lowpass  filter  for  the  STFT  is  set  to  give  unit  power  throughput,  the  noise  power 
at  the  output  of  the  filterbank  is  identical  to  the  spectral  estimates  described  above. 
Since  the  noise  is  assumed  to  be  spatially  white,  the  noise  power  in  the  estimate  of 
mode  m  is  simply  the  estimated  noise  level  for  that  frequency  bin  times  the  squared 
length  of  the  weight  vector  for  that  mode,  i.e.,  w“[cu*]wm[<Uifc]. 

In  the  numerical  results  presented  in  this  chapter,  the  plotting  thresholds  are  set 
on  a  reception-by-reception  basis  as  follows.  For  plots  of  the  received  pressure  time 
series,  the  0  dB  power  level  corresponds  to  the  estimated  noise  floor  associated  with 
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the  75  Hz  bin.  In  the  modal  time  series,  the  0  dB  level  corresponds  to  the  noise  level 
in  mode  1  at  75  Hz  for  that  reception. 

4.2  Simulated  Data  Set 

For  comparison  purposes,  a  simulated  data  set  has  been  generated  using  the  parabolic 
equation  code  RAM  [14].  The  simulations  are  designed  to  model  propagation  along 
the  Pioneer-Hawaii  path.  Background  sound  speed  profiles  for  the  simulation  en¬ 
vironment  are  the  Levitus-winter  profiles  along  the  geodesic  path  between  Pioneer 
Seamount  and  the  Hawaii  array,  shown  in  Chapter  2.  The  bathymetry  for  the  sim¬ 
ulations  corresponds  to  actual  bathymetry  along  the  path  with  the  exception  of  one 
simplification:  the  steep  downslope  near  the  source  has  been  eliminated.  The  source 
is  at  939.5  meters  depth  (the  depth  of  the  seamount  where  the  source  is  located); 
bottom  depth  at  the  source  range  is  3317  m.  Section  4.3.3  illustrates  the  effects  of 
adding  the  seamount  back  into  the  simulations,  and  indicates  that  source  bathymetry 
may  affect  the  spread  of  the  final  arrivals. 

The  simulated  data  set  consists  of  multiple  realizations  of  the  received  time  series 
on  a  40-element  array  with  sensors  located  at  the  nominal  ATOC  VLA  depths.  In 
addition  to  one  time  series  for  the  Levitus  background  environment,  there  are  15 
time  series  which  include  sound  speed  perturbations  due  to  internal  waves.  Real¬ 
izations  of  the  internal  wave  (IW)  field  were  generated  using  the  method  of  Colosi 
and  Brown  [13].  A  1/2  Garrett-Munk  strength  was  used  for  all  simulations,  in  ac¬ 
cordance  with  the  earlier  findings  of  Colosi.  Ten  of  the  15  IW  simulations  are  for 
independent  realizations  of  the  internal  wave  field.  The  other  5  are  for  successive 
time  steps  associated  with  one  of  the  independent  realizations. 

The  PE  code  generates  solutions  of  the  frequency-dependent  wave  equation. 
Broadband  simulations  therefore  require  Fourier  synthesis  of  the  PE  results  to  pro¬ 
duce  a  time  series.  The  frequency-spacing  used  for  the  PE  calculations  determines 
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the  time  interval  that  can  be  represented  without  aliasing.  For  the  ATOC  environ¬ 
ment,  the  time  spread  of  all  the  arrivals  (ray-like  and  mode-like)  is  on  the  order  of 
10-20  seconds,  implying  a  required  frequency  spacing  of  0.1-0.05  Hz.  Since  the  low 
mode  arrivals  are  the  signals  of  primary  interest  in  this  thesis,  a  simplification  is 
possible,  which  permits  a  wider  frequency  spacing  (thus  a  significant  computational 
savings).  In  these  simulations  the  pressure  field  is  calculated  for  each  frequency  and 
projected  onto  the  modes  of  the  receiver  environment.  Then,  only  a  subset  of  those 
modes  are  used  to  synthesize  the  time  series.  Empirical  results  indicate  that,  for  this 
frequency  range,  45  modes  are  fully  contained  within  an  8  second  interval;  25  modes 
are  contained  within  a  4  second  window.  Thus,  a  frequency  spacing  of  0.125  Hz  or 
0.25  Hz  is  permitted  provided  that  the  time  series  is  synthesized  with  only  45  or 
25  modes,  respectively.  Note  that  4-second  simulations  are  more  than  adequate  for 
examining  the  behavior  of  the  first  10  modes. 

The  simulations  do  not  contain  any  additive  noise  components.  In  order  to  use 
a  dB  scaling  for  the  plots  that  is  consistent  with  the  real  data,  it  is  assumed  that 
there  is  a  noise  floor  12  dB  below  the  peak  in  the  pressure  field.  This  level  was 
chosen  based  on  the  average  peak-to-noise  ratio  in  Figure  4-4.  Simulated  pressure 
time  series  are  scaled  so  that  this  assumed  noise  floor  is  at  0  dB;  mode  time  series 
are  scaled  with  respect  to  the  assumed  noise  floor  for  the  mode  1  estimate. 

4.3  Modal  Time  Series:  Experiment  vs. 
Simulation 

This  section  examines  the  mode  estimates  for  one  ATOC  reception  and  compares 
them  to  estimates  for  two  simulated  receptions.  One  simulation  is  for  the  Levitus 
background  environment;  the  other  is  for  the  background  environment  plus  sound 
speed  perturbations  due  to  internal  waves.  The  purpose  of  this  section  is  to  high¬ 
light  important  features  of  the  short-time  mode  spectra  and  to  illustrate  qualitative 
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agreement  between  the  ATOC  data  and  the  simulation  that  includes  internal  wave 
effects. 

4.3.1  ATOC  Reception 

Figure  4-5  is  a  plot  of  a  demodulated  pressure  time  series  (one  4-period  average) 
recorded  on  the  Hawaii  array  in  late  December  1995.  For  reference,  the  pressure  time 
series  for  the  simulations  discussed  later  are  displayed  below  the  ATOC  reception 
(Figures  4-6  and  4-7).  In  the  ATOC  reception,  the  early-arriving  planewaves  (or 
ray  arrivals)  are  rather  weak.  Unlike  the  adiabatic  examples  considered  in  previous 
chapters,  there  are  no  immediately  identifiable  modes  contained  in  the  late-arriving 
energy. 

The  mode  estimates  for  this  reception  are  calculated  using  short-time  Fourier 
techniques  developed  in  Chapter  3.  Recall  that  the  STFT  processor  designed  for  the 
ATOC  data  uses  an  0.4  second  Hanning  window  and  a  10-mode  pseudo-inverse  spatial 
filter.  Note  that  for  the  plots  in  this  chapter,  the  STFT  estimates  are  computed  for 
the  60-90  Hz  band,  with  a  bin  spacing  of  1.25  Hz.  The  latter  is  chosen  for  plotting 
purposes  only;  the  critical  frequency  spacing  required  by  the  choice  of  window  is 
2.5  Hz. 

Figure  4-8  shows  the  frequency-stacked  plots  for  the  first  10  modes  of  this  ATOC 
reception.  There  are  a  number  of  observations  to  make  about  these  results.  First, 
each  mode  consists  of  multiple  arrivals,  spread  over  2-3  seconds.  This  is  a  stark 
contrast  to  the  single,  dispersive  arrivals  that  characterized  the  adiabatic  propagation 
example  discussed  in  the  last  chapter.  Also,  recall  that  in  the  adiabatic  case,  the 
modes  arrive  in  descending  order  and  the  highest  modes  are  temporally  separable 
from  the  lowest  modes.  For  this  ATOC  reception,  there  is  no  obvious  ordering  of  the 
arrivals  and  the  spread  of  the  signals  is  such  that  the  first  ten  modes  overlap  in  time. 

The  third  point  to  consider  concerns  a  common  feature  of  multipath  channels: 
frequency-selective  fading.  This  type  of  fading  occurs  when  two  signals  with  different 
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Figure  4-5:  ATOC  Reception 


Figure  4-6:  PE  simulation  without  internal  waves 
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Figure  4-7:  PE  simulation  with  internal  waves  at  1/2  Garrett-Munk  strength 
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phase  characteristics  arrive  at  the  same  time.  The  spectrum  for  the  superposition 
of  these  signals  may  contain  amplitude  fades  due  to  destructive  interference.  If  the 
individual  arrivals  are  temporally-resolvable  using  a  short-time  transform,  then  the 
characteristics  for  each  path  can  be  measured.  In  this  case,  the  magnitude  across 
frequency  should  not  contain  deep  fades. 

The  STFT  processor  used  on  the  ATOC  data  can  begin  to  temporally  resolve 
signals  at  a  separation  of  0.2  seconds  (see  discussion  on  page  69).  Note  that  since 
bins  with  2.5  Hz  spacing  are  slightly  correlated  (due  to  the  overlap  of  the  bandpass 
filters  associated  with  the  STFT),  signals  must  show  coherence  over  bandwidths 
of  greater  than  5  Hz  before  it  is  clear  that  the  coherence  is  due  to  the  underlying 
arrival  and  not  an  artifact  of  processing.  Although  the  plots  in  Fig.  4-8  indicate 
the  presence  of  faded  arrivals,  they  also  contain  some  arrivals  that  extend  over  5- 
10  Hz  bandwidths.  This  suggests  that  the  STFT  processor  may  be  resolving  some 
individual  multipaths. 

4.3.2  Simulated  Receptions 

The  previous  section  reviewed  important  characteristics  of  the  mode  arrival  structure 
estimated  from  ATOC  data.  For  comparison,  this  section  considers  two  of  the  simu¬ 
lated  receptions  (described  in  Section  4.2).  The  first  is  for  propagation  through  the 
Levitus  background  environment;  the  second  includes  internal-wave-induced  sound 
speed  perturbations. 

Figure  4-6  shows  the  received  time  series  for  propagation  through  the  background 
environment.  Note  that  individual  mode  arrivals  are  evident  in  the  pressure  field  e.g., 
mode  2  is  the  strongest  final  arrival  located  at  2375  seconds.  (Mode  1  is  not  strongly 
excited  by  a  source  located  at  939.5  m  depth.)  Figure  4-9  displays  the  short-time 
Fourier  mode  estimates  for  the  first  10  modes,  indicating  that  these  modes  arrive 
in  descending  order.  The  solid  black  line  on  each  of  the  frequency-stacked  plots 
represents  the  adiabatic  arrival  times.  Agreement  between  the  adiabatic  predictions 
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and  the  PE  simulation  without  internal  waves  is  quite  good.  While  this  does  not 
prove  that  propagation  through  the  Levitus  background  environment  is  adiabatic,  it 
does  demonstrate  that,  in  the  absence  of  internal  waves,  each  of  the  first  10  modes 
is  dominated  by  a  single  dispersive  arrival. 

In  contrast  to  the  background  simulation,  the  results  for  the  internal  wave  simu¬ 
lation  show  more  complicated  mode  arrival  patterns.  Figure  4-7  is  the  pressure  time 
series  for  the  internal  wave  simulation,  and  Figure  4-10  displays  the  corresponding 
mode  estimates.  Similar  to  the  ATOC  data,  the  estimates  for  the  perturbed  envi¬ 
ronment  contain  multiple  arrivals  in  each  mode.  Some  of  these  arrivals  also  exhibit 
frequency-selective  fading,  and  others  appear  to  be  coherent  over  bands  on  the  order 
of  10  Hz. 

To  permit  a  comparison  between  experiment  and  simulation,  Figure  4-11  shows 
the  estimated  spectra  for  modes  1  and  10  for  the  ATOC  reception  and  the  two  sim¬ 
ulated  receptions.  For  reference,  adiabatic  travel  time  predictions  for  the  California- 
Hawaii  path  are  shown  at  the  bottom.  There  is  qualitative  similarity  between  the 
mode  estimates  for  the  ATOC  reception  and  the  PE  simulation  containing  internal 
waves,  i.  e.,  both  contain  multipath  arrivals  and  show  evidence  of  frequency-selective 
fading.  Although  mode  10  starts  arriving  before  mode  1  in  both  the  ATOC  data 
and  the  internal  wave  simulation,  there  is  significant  overlap  of  these  two  modes  in 
time.  This  is  in  contrast  to  the  background  simulation  and  the  adiabatic  case  where 
those  two  modes  are  temporally  separated.  Note  that  the  ATOC  data  does  show 
significantly  more  time-spread  than  the  simulation. 

The  most  striking  disagreement  between  the  PE  simulations  and  the  real  data 
is  that  the  simulations  have  a  sharp  cutoff  around  2375  seconds  (obvious  in  both 
the  pressure  time  series  and  the  mode  estimates),  but  there  is  no  discernible  cutoff 
in  the  ATOC  receptions.  In  the  ATOC  data,  high  energy  arrivals  occur  before 
2375.5  seconds,  but  the  signal  trails  off  with  a  series  of  lower  energy  peaks  after  that. 
Heaney  has  associated  this  “afterglow”  in  the  arrival  pattern  with  bottom  interaction 
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Figure  4-11:  Comparison  of  ATOC  data,  PE  simulations,  and  adiabatic  predictions 
for  modes  1  and  10. 


near  the  source  [53].  The  next  section  considers  the  downslope  propagation  issue 
briefly. 

4.3.3  Downslope  Propagation  Example 

Although  a  thorough  analysis  of  the  downslope  propagation  problem  is  beyond  the 
scope  of  this  thesis,  it  is  useful  to  consider  one  example  which  illustrates  the  effects 
of  the  slope  on  the  initial  excitation  of  the  modes. 

Figure  4-12  shows  the  actual  bathymetry  (dash-dot  line)  measured  during  source 
deployment,  along  with  several  approximations.  The  dotted  line  represents  the  bot¬ 
tom  depths  used  in  the  simulations  described  in  Section  4.2.  Those  ignore  the  slope 
entirely.  As  indicated  by  the  figure,  the  seamount  has  very  steep  sides;  at  its  steepest, 
the  slope  off  Pioneer  is  approximately  20  degrees. 


First  consider  a  PE  simulation  with  the  actual  sloping  bottom  near  the  source. 
The  PE  code  does  not  handle  shear  waves,  but  a  reasonable  approximation  (ignoring 
interface  waves)  for  a  solid  bottom  is  to  use  a  fluid  bottom  with  a  compressional 
speed  equal  to  that  of  the  shear  parameters  in  the  actual  bottom  [4],  For  the  ATOC 
source  site,  the  bottom  is  basalt  which  has  a  shear  speed  of  2500  m/s,  a  density 
of  2.7  g/cm3,  and  an  attenuation  of  0.1  dB  per  wavelength.  These  are  the  bottom 
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parameters  used  in  all  of  the  simulations.  Figure  4-13  compares  the  received  pressure 
time  series  for  the  actual  slope  environment  (top  plot)  and  the  no  slope  environment 
(bottom  plot).  These  two  simulations  include  the  same  internal  wave  perturbations. 
The  plots  are  scaled  to  the  peak  pressure  over  the  two  receptions.  Note  that  the 
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Figure  4-13:  Comparison  of  propagation  through  the  environment  that  includes  the 
actual  sloping  bottom  near  the  source  (top  plot)  and  the  environment  without  the 
slope  (bottom  plot) 

presence  of  a  slope  changes  the  final  arrivals  and  adds  approximately  1  second  of 
additional,  lower  amplitude  arrivals.  Figure  4-14  illustrates  the  effects  of  the  slope 
on  the  mode  1  and  mode  10  estimates.  Additional  arrivals  are  added  to  both  modes 
after  2375  seconds.  There  is  also  a  strong  arrival  added  in  mode  10,  before  the  cutoff. 
Although  mode  1  does  not  appear  to  be  as  strongly  affected  as  mode  10  in  this  case, 
it  is  unclear  whether  that  would  be  true  for  other  realizations  of  the  internal  wave 
field. 

To  get  insights  into  what  is  happening  near  the  source,  consider  a  plot  of  the 
time  series  for  the  first  80  modes  at  a  distance  of  50  km  from  the  array.  Figure  4-15 


Time  (seconds)  Time  (seconds) 


Figure  4-14.  Comparison  of  arrivals  in  modes  1  and  10  for  the  sloping  bottom  (top 
plots)  and  the  zero-slope  approximation  (bottom  plots) 

shows  these  results  for  the  actual  slope  and  for  the  3  approximations  to  that  slope 
in  order  to  illustrate  how  steep  the  slope  has  to  be  to  produce  this  effect.  These 
are  broadband  time  series  computed  by  projecting  onto  the  modes  at  50  km  at  each 
frequency  and  then  synthesizing  a  time  series  for  each  mode.  The  top  plot  shows  the 
results  for  a  simulation  with  the  actual  source  bathymetry.  In  most  of  the  modes 
below  60,  the  energy  is  concentrated  in  two  pulses  that  are  spaced  approximately 
one  second  apart.  This  is  in  contrast  to  the  no  slope  case,  shown  in  the  bottom 
plot,  which  contains  only  one  dominant  pulse  in  each  mode.  The  two  middle  plots 
indicate  how  the  second  pulse  fades  out  as  the  slope  is  decreased.  The  presence 
of  the  seamount  clearly  changes  the  mode  excitation,  but  more  work  is  needed  to 
verify  the  actual  mode  scattering  mechanisms  down  the  slope  that  produce  this 
effect.  Heaney  attributes  the  afterglow  to  a  bottom  bounce  near  the  source  [53], 
It  seems  likely  that  the  second  pulse  in  the  mode  excitation  is  due  to  that  bottom 
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Figure  4-15:  Comparison  of  the  time  series  for  the  first  80  modes  at  50  km  range 
from  the  Pioneer  Seamount  source.  Top  plot  is  for  the  actual  slope.  Middle  two  plots 
are  for  the  8.5  degree  and  4.3  degree  approximations,  respectively.  The  bottom  plot 
corresponds  to  the  no  slope  case. 
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bounce.  Although  more  research  is  clearly  required  in  this  area,  these  simulation 
results  indicate  that  the  ATOC  measurements  may  represent  the  response  of  the 
ocean  channel  to  a  superposition  of  two  temporally-separated  pulses  rather  than  the 
response  to  a  single  pulse. 

This  section  (4.3)  has  examined  mode  estimates  for  a  single  ATOC  reception  and 
compared  them  to  results  for  a  simulation  through  internal  waves  at  1/2  Garrett- 
Munk  strength.  The  ATOC  results  and  the  simulation  agree  qualitatively,  with 
both  indicating  that  the  mode  signals  at  megameter  ranges  contain  multiple,  faded 
arrivals  rather  than  the  single  dispersive  arrivals  that  characterize  adiabatic  propaga¬ 
tion.  Although  the  similarities  between  experiment  and  simulation  are  encouraging, 
the  ATOC  data  does  show  significantly  more  time-spread.  The  most  striking  dis¬ 
agreement  between  the  simulated  time  series  and  the  measured  data  is  the  sharp 
cutoff  in  the  simulation,  which  was  shown  to  be  a  consequence  of  ignoring  the  source 
bathymetry. 

4.4  Temporal  Variability 

The  purpose  of  this  section  is  to  consider  the  temporal  variability  of  the  mode  arrivals. 
In  previous  work,  Tang  and  Tappert  obtained  measurements  of  the  temporal  coher¬ 
ence  of  pulses  propagating  through  internal  waves  in  a  shallow  (200  m  deep)  water 
environment  at  ranges  of  20  km.  They  considered  the  coherence  times  of  pulses 
and  did  not  address  the  issue  of  modes  specifically.  For  deep  water  environments, 
Colosi  et  al.  [71]  (citing  the  work  of  Flatte  and  Stoughton  [72])  indicate  that  acoustic 
coherence  times  are  on  the  order  of  tens  of  minutes,  whereas  the  coherence  times  for 
internal  waves  are  on  the  order  of  hours.  In  analyzing  data  from  the  ATOC  Engi¬ 
neering  Test  (range=3250  km),  Worcester  et  al.  used  12.7  minute  coherent  averages 
for  the  ray  arrivals,  since  the  SNR  did  not  increase  for  longer  averaging  times.  For 
that  same  data  set,  Colosi  et  al.  concluded  that  time  fluctuations  in  the  wavefronts 
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show  no  coherence  at  2  hour  lag  times.  Note  that  in  the  works  cited  above,  the  focus 
was  primarily  on  analyzing  the  ray  arrivals,  rather  than  the  late-arriving  modes.  The 
rest  of  this  section  considers  the  temporal  variations  of  the  mode  arrivals  in  ATOC. 

Recall  that  for  each  source  transmission,  the  ATOC  VLA  recorded  10  four-period 
averages  of  the  27.28  second  pseudo-random  sequence.  The  results  discussed  in 
Section  4.3.1  were  computed  using  the  first  4-period  average  of  reception  363004024 
at  the  Hawaii  array.  Figure  4-16  compares  mode  estimates  for  the  first  4-period 
average  and  the  last  4-period  average  for  that  reception.  The  top  plots  are  the 
estimated  spectra  for  modes  1  and  10  of  the  first  average  in  the  transmission  and  the 
bottom  plots  correspond  to  the  last  average.  Note  that  the  time  difference  between 
the  end  of  the  first  reception  and  the  start  of  the  tenth  is  approximately  14.5  minutes. 
The  plot  clearly  demonstrates  that  the  mode  signals  change  significantly  over  that 
time  interval:  some  arrivals  drop  out  entirely,  and  new  arrivals  emerge. 
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Figure  4-16:  Comparison  of  modes  1  and  10  for  the  first  (top  plots)  and  last  (bottom 
plots)  periods  of  a  source  transmission 


For  comparison,  consider  the  results  of  PE  simulations  that  take  into  account 
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Figure  4-17:  Comparison  of  modes  1  and  10  for  realizations  of  a  time- varying  internal 
wave  environment.  The  lag  time  between  the  top  and  bottom  plots  is  16  minutes. 


time  variations  of  the  internal  wave  field.  Using  Colosi  and  Brown’s  approach,  a 
realization  of  the  time-varying  internal  wave  field  was  generated  for  a  nominal  time 
to,  and  for  to+16  minutes.  Fig.  4-17  shows  the  estimates  for  modes  1  and  10  that 
result  from  these  simulations.  These  plots  demonstrate  fade-ins  and  fade-outs  that 
are  similar  to  those  seen  in  the  ATOC  data. 

Based  on  the  single  realizations  of  the  ATOC  data  and  the  simulated  data  shown 
above,  the  low  order  mode  arrivals  fluctuate  over  time  intervals  on  the  order  of  min¬ 
utes.  The  following  approach  was  used  to  obtain  a  measure  of  the  average  coherence 
time  for  each  mode.  First,  correlation  coefficients  were  computed  at  a  set  of  10  lags 
for  the  received  signals  in  a  particular  mode,  corresponding  to  the  10  four-period 
averages.  Note  that  the  time  interval  used  for  these  calculations  was  2373  seconds 
to  2376  seconds,  where  the  mode  arrivals  are  concentrated.  These  correlation  coef¬ 
ficients  were  then  averaged  over  96  ATOC  receptions  (from  yearday  363  to  yearday 
435).  The  results  are  shown  in  Figure  4-18.  The  left  plot  in  the  figure  shows  the 
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Figure  4-18:  Temporal  coherence  as  a  function  of  frequency  for  mode  1  (left  plot) 
and  as  a  function  of  mode  number  for  the  75  Hz  bin  (right  plot).  These  results  were 
obtained  by  averaging  across  96  receptions  (yeardays  363-435). 

temporal  coherence  of  mode  1  for  three  different  frequencies.  There  does  appear 
to  be  some  mild  frequency  dependence,  with  mode  1  in  the  65  Hz  bin  being  more 
correlated  than  the  same  mode  at  85  Hz.  For  mode  1,  correlation  decreases  to  a  level 
of  0.5  between  6  and  8  minutes.  The  right  plot  in  the  figure  shows  coherence  as  a 
function  of  mode  number  for  the  75  Hz  bin.  Based  on  these  results,  temporal  coher¬ 
ence  is  not  a  function  of  mode  number,  which  is  not  surprising  given  the  amount  of 
mode  coupling. 

The  above  calculations  provide  a  measure  of  the  coherence  time  of  the  whole 
signal  in  a  particular  mode,  z.e.,  all  the  multipath  arrivals.  An  interesting  question 
to  ask  is  whether  some  peaks  are  more  stable  than  others  across  the  18.2  minute 
transmission  interval.  Figure  4-19  illustrates  how  mode  signals  in  the  75  Hz  bin  vary 
over  the  18.2  minutes.  The  plots  are  stacks  of  the  received  signals  in  mode  1  (top)  and 
mode  10  (bottom)  at  successive  1.82  minute  lags.  Some  of  the  peaks  are  consistent 
across  the  full  transmission  interval,  while  others  fade  out  suddenly,  e.g.,  the  mode  1 
arrival  at  2375  seconds  that  drops  out  around  12  minutes.  Again,  it  is  interesting  to 
compare  these  results  to  PE  simulations  involving  time-varying  internal  wave  fields. 
For  these  results,  the  internal  wave  sound  speed  perturbations  were  computed  at 
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4-minute  intervals  between  0  and  20  minutes.  Figure  4-20  shows  the  resulting  stacks 
for  modes  1  and  6  for  the  75  Hz  bin.  This  data  confirms  that  the  PE  simulations 
with  internal  waves  at  1/2  Garrett  Munk  strength  result  in  qualitatively  similar  time 
variations  of  the  modes. 
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Figure  4-19:  Variability  of  modes  1  and  10  across  a  single  transmission  (18.2  minutes) 
in  the  75  Hz  bin 


This  section  has  shown  that  the  mode  arrival  structure  varies  considerably  over 
intervals  on  the  order  of  20  minutes.  The  following  section  considers  four  hour 
intervals  and  asks  whether  averaging  over  multiple  transmissions  can  be  used  to 
obtain  insight  into  mode  behavior. 


4.5  Mode  Statistics 


On  days  that  the  ATOC  source  was  on  there  were  transmissions  at  4-hour  intervals. 
Given  the  large  amount  of  variability  in  the  mode  signals  over  a  20-minute  period,  it 
is  expected  that  the  transmissions  at  four-hour  intervals  will  sample  very  different  re- 
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Figure  4-20:  Variability  of  modes  1  and  6  for  the  simulated  data  in  the  75  Hz  bin, 
calculated  at  4  minute  intervals. 

alizations  of  the  internal  wave  field.  Such  behavior  can  be  seen  in  Figure  4-21,  which 
shows  the  short-time  spectra  for  mode  1  (the  first  4-period  average)  for  three  consec¬ 
utive  transmissions.  Clearly,  the  mode  arrival  structure  is  significantly  different  at 
4-hour  lags.  Although  all  arrivals  for  all  three  receptions  shown  are  approximately 
centered  around  2375  seconds,  and  there  is  some  consistency  in  the  weaker  peaks 
( e.g .,  the  faded  arrival  around  2376  seconds),  the  strong  arrivals  do  not  occur  in  the 
same  places  for  each  reception,  and  overall  the  interference  pattern  is  quite  different. 

To  get  an  idea  of  how  the  peak  arrivals  are  distributed  as  a  function  of  time, 
consider  the  histograms  shown  in  Figure  4-22.  These  plots  were  compiled  from  the 
first  96  good  receptions  recorded  on  yeardays  363  through  436,  using  a  detection 
threshold  of  12  dB  above  the  estimated  noise  level  for  each  mode.  These  results 
confirm  that  there  is  not  a  dominant  arrival  time  in  each  mode;  rather,  the  peaks 
are  distributed  between  2373  seconds  and  2376  seconds  for  all  of  the  modes.  Note 
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Figure  4-21:  Mode  1  in  ATOC  receptions  at  4-hour  intervals 


that  the  distribution  of  the  higher  modes  ( e.g .,  mode  10)  is  skewed  towards  the  early 
part  of  the  3-second  interval  while  the  low  modes  are  concentrated  towards  the  latter 
part  of  the  interval.  This  is  consistent  with  the  deep  water  dispersion,  where  higher 
modes  arrive  first. 

Given  the  amount  of  variability  in  the  arrival  structure,  does  averaging  over 
multiple  transmissions  reveal  consistent  features  in  the  mode  signals?  To  answer 
this  question,  it  is  useful  to  consider  several  statistics:  leading  edge,  falling  edge,  and 
centroid.  Note  that  the  leading  and  falling  edges  are  (respectively)  the  first  and  last 
time  indices  where  the  complex  envelope  of  the  mode  estimate  exceeds  a  threshold. 
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Figure  4-22:  Histogram  of  peak  arrivals  in  the  75  Hz  bin.  Threshold  for  peak  detec¬ 
tion  was  set  at  12  dB  above  the  noise  floor.  These  results  were  computed  from  the 
first  96  good  receptions  (yeardays  363  to  436). 


The  centroid  is  defined  to  be  the  center  of  mass  of  the  portion  of  the  complex  envelope 
that  is  higher  than  the  threshold  value.  For  all  of  the  results  presented  below,  the 
threshold  was  set  at  12  dB  above  the  noise  floor  in  each  mode.  In  order  to  examine 
how  these  statistics  vary  as  a  function  of  time  over  the  course  of  the  experiment,  the 
ATOC  receptions  were  divided  into  13  groups  as  indicated  in  Fig.  4-1. 

Figure  4-23  shows  the  leading  edge,  falling  edge,  and  centroid  in  the  75  Hz  bin  for 
the  first  group  of  receptions.  This  group  consists  of  20  receptions,  each  containing 
10  four-period  averages.  The  data  is  plotted  so  that  each  whole  number  on  the  x- 
axis  represents  the  start  of  a  set  of  4-period  averages.  As  this  figure  indicates,  the 
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Figure  4-23:  Leading  edge,  falling  edge,  and  centroid  in  the  75  Hz  bin  for  the  first 
20  ATOC  receptions,  each  consisting  of  10  four-period  averages. 

three  statistics  of  interest  can  fluctuate  significantly  over  one  source  transmission 
(18.2  minutes)  as  arrivals  fade  in  and  out. 

Figure  4-24  shows  the  results  of  averaging  over  all  the  receptions  in  the  first 
group  (a  total  of  2003)  to  obtain  the  leading  and  falling  edges  as  a  function  of 
frequency  for  the  first  10  modes.  For  reference,  the  figure  also  includes  the  average 
leading  and  falling  edges  for  the  mode  estimates  of  10  simulated  receptions.  The 

3Note  that  the  average  may  contain  less  than  200  points  since  there  are  a  few  receptions  where 
no  part  of  the  estimate  exceeds  the  detection  threshold. 
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ATOC  data:  average  leading  and  falling  edges  for  group  1 


Mode  1 
Mode  2 
Mode  3 
Mode  4 
Mode  5 
Mode  6 
Mode  7 
Mode  8 
Mode  9 
Mode  1 0 


2373.5 


2374  2374.5  2375 

Time  (seconds) 

Simulation  data:  average  leading  and  falling  edges 
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Figure  4-24:  Comparison  of  average  leading  and  falling  edges  for  the  first  group  of 
ATOC  receptions  and  the  simulated  data  set 
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error  bars  on  these  plots  (and  all  subsequent  ones)  represent  the  estimated  RMS 
error,  he.,  the  sample  standard  deviation  divided  by  the  square  root  of  the  number 
of  samples  included  in  the  average.  The  ATOC  statistics  reveal  several  significant 
features.  First,  the  average  spread  between  leading  and  falling  edges  is  on  the  order  of 
1.5  seconds.  Second,  the  high  modes  have  earlier  arrival  times,  while  low  modes  have 
later  arrival  times,  as  would  be  expected  in  a  deep  water  channel.  Although  there  is 
some  crossover  between  neighboring  modes,  the  leading  and  falling  edges  show  that 
there  are  statistically  significant  differences  in  arrival  time  among  the  modes,  e.g ., 
compare  modes  1  and  10.  Finally,  note  that  the  plot  indicates  that  the  signals  at  the 
center  frequency  (75  Hz)  are  more  spread  than  those  at  either  end  of  the  band.  This 
effect  may  be  partially  due  to  the  fact  that  the  same  detection  threshold  is  used  for 
all  frequency  bins,  while  the  source  spectrum  rolls  off  as  a  function  of  frequency. 

Since  the  simulated  data  is  averaged  over  only  ten  receptions,  the  resulting  curves 
in  Fig.  4-24  are  not  as  smooth  as  the  ATOC  data  and  the  error  bars  are  larger. 
Nevertheless,  the  simulated  data  is  comparable  to  the  real  data  in  two  important 
respects:  the  leading  edges  show  similar  behavior  (as  a  function  of  frequency)  to  the 
ATOC  data,  and  the  arrival  times  are  also  comparable.  Unlike  the  falling  edges  in 
the  ATOC  data,  however,  the  falling  edges  in  the  simulated  data  are  much  more 
concentrated.  This  sudden  cutoff  has  been  noted  previously  and  is  likely  the  result 
of  the  lack  of  modeling  of  the  downslope  propagation  in  the  simulation. 

Now  consider  the  average  centroid  statistics  as  a  function  of  frequency.  Figure  4- 
25  compares  the  average  centroid  locations  for  the  first  group  of  ATOC  receptions 
with  the  average  centroids  for  the  simulated  receptions.  Predicted  adiabatic  travel 
times  are  included  for  reference.  The  experimental  data  centroid  data  shows  the 
modes  arriving  in  descending  order  (though  there  is  some  crossover  among  nearest- 
neighbors).  For  the  highest  modes,  there  is  a  slight  trend  across  frequency:  the  higher 
frequencies  arrive  later,  as  would  be  expected  in  deep  water.  The  simulated  data  show 
good  agreement:  the  centroids  extend  over  approximately  the  same  time  interval 
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Figure  4-25:  Comparison  of  average  centroid  locations  for  the  first  group  of  ATOC 
receptions,  the  simulated  receptions,  and  adiabatic  predictions 
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(~  0-8  seconds),  and  there  is  again  a  slight  trend  of  increasing  arrival  time  with 
increasing  frequency  in  the  higher  modes.  This  supports  the  idea  that  1/2  Garrett- 
Munk  strength  is  appropriate.  The  key  difference  between  the  ATOC  data  and  the 
simulated  data  in  Figure  4-25  is  a  shift  in  the  mean  arrival  times:  the  centroids 
for  the  first  group  of  ATOC  data  occur  approximately  0.2  seconds  later  than  the 
centroids  of  the  simulated  receptions.  Comparing  the  ATOC  and  simulated  data 
to  the  predicted  adiabatic  arrival  times,  one  key  difference  is  evident:  the  interval 
over  which  the  arrival  times  are  spread  is  much  longer  in  the  adiabatic  case.  This  is 
expected  because  mode  scattering  accounts  for  the  concentration  in  the  ATOC  and 
simulated  data. 

Figures  4-26  and  4-27  show  the  centroids  as  a  function  of  frequency  for  the  13 
groups  of  receptions  defined  in  Figure  4-1.  These  receptions  were  taken  over  a  period 
of  approximately  5  months,  between  the  end  of  December  of  1995  and  May  of  1996. 
The  centroids  in  each  of  the  groups  show  similar  behavior  as  a  function  of  frequency, 
but  there  is  clearly  a  drift  in  the  arrival  times  from  one  set  of  receptions  to  the 
next.  The  drift  in  arrival  times  is  more  obvious  in  the  plot  of  Figure  4-28,  which 
shows  the  average  centroids  of  first  10  modes  at  75  Hz  as  a  function  of  yearday. 
Note  that  the  minimum  travel  time  for  all  of  the  modes  occurs  around  yearday  427, 
which  corresponds  to  the  beginning  of  March,  1996.  The  trend  of  decreasing  mode 
arrival  time  from  the  start  of  the  experiment  until  March  and  increasing  afterwards 
agrees  with  the  trend  observed  in  the  ray  arrivals  for  the  Hawaii  VLA  [73].  The 
difference  in  the  arrival  times  between  the  first  group  of  receptions  (yearday  363) 
and  the  sixth  group  of  receptions  (yearday  427)  is  between  0.3  and  0.5  seconds  for 
the  first  10  modes.  A  T-test  confirms  that  these  time  differences  are  statistically 
significant  (using  a  significance  level  of  0.01).  Based  on  the  1998  Science  article  by 
the  ATOC  Consortium  [73],  the  ray  arrivals  exhibit  somewhat  smaller  time  shifts  on 
the  order  of  0.25  seconds  over  that  same  period.  The  Science  article  notes  that  the 
trend  observed  at  the  Hawaii  array  is  not  in  agreement  with  the  expected  seasonal 
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Group  1 


Group  5 


Time  (seconds) 

Figure  4-26:  Average  centroids  as  a  function  of  frequency  for  the  first  7  groups  of 
ATOC  receptions.  Legend  is  identical  to  that  shown  in  Fig.  4-25.  See  Fig.  4-1  for  a 
definition  of  the  groups. 
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Figure  4-27:  Average  centroids  as  a  function  of  frequency  for  groups  8-13  of  the 
ATOC  receptions.  Legend  is  identical  to  that  shown  in  Fig.  4-25.  See  Fig.  4-1  for  a 
definition  of  the  groups. 
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Figure  4-28:  Average  centroids  in  the  75  Hz  bin  for  modes  1-10  as  a  function  of 
yearday 


trend  (which  would  have  shown  increasing  travel  times  in  winter  and  decreasing  in 
the  summer)  and  postulates  that  this  is  due  to  a  subsurface  warming  near  the  receiver 
that  offsets  the  winter  surface  cooling  layer  near  the  source. 

The  trend  in  arrival  times  as  a  function  of  yearday  is  apparent  in  the  falling 
edges  of  the  modes,  but  is  not  as  evident  in  the  leading  edge  statistic.  Figure  4-29 
illustrates  this  point  using  mode  1  in  the  75  Hz  bin  as  an  example.  These  curves 
are  representative  of  the  results  for  other  low  modes.  Note  that  the  falling  edge 
shows  a  decrease  in  arrival  time  from  January  to  March  that  correlates  well  with  the 
centroid  data.  The  maximum  shift  in  the  mean  falling  edge  time  is  also  0.45  seconds. 
While  the  leading  edge  data  show  a  significant  decrease  in  arrival  time  for  the  set  of 
receptions  around  yearday  427,  there  is  no  obvious  downward  trend  in  the  first  five 
groups  of  receptions  as  there  is  in  the  centroid  and  falling  edge  statistics. 

It  is  important  to  consider  the  errors  associated  with  these  statistics.  For  the 
centroids  and  falling  edges,  Figure  4-30  shows  the  average  (over  the  13  groups  of 
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receptions)  RMS  error  at  75  Hz  for  the  first  10  modes.  The  average  RMS  error 
for  the  centroids  is  between  20-27  milliseconds,  with  a  slight  trend  of  increasing 
error  with  mode  number.  Error  for  the  falling  edge  statistics  tends  to  be  larger 
(29-34  milliseconds)  and  does  not  exhibit  a  mode-dependent  trend.  Note  that  the 


RMS  error  averaged  over  13  groups  of  ATOC  receptions 


75  Hz 


Figure  4-30:  Average  of  the  RMS  error  over  13  groups  of  receptions  as  a  function  of 
mode  number.  Results  are  shown  for  the  centroids  and  the  falling  edges  of  the  75  Hz 
bin. 

order  of  magnitude  of  these  errors  is  comparable  to  the  measured  fluctuations  of  the 
identifiable  ray  arrivals  (11-19  ms)  and  the  pulse  termination  (22  ms)  [71]. 


4.6  Summary 

This  chapter  has  provided  the  first  detailed  look  at  the  low-mode  signals  in  the 
ATOC  receptions  at  Hawaii.  While  much  work  remains  to  be  done,  some  impor¬ 
tant  conclusions  about  the  nature  of  the  low-mode  arrivals  can  be  drawn  from  this 
analysis. 

The  first  observation  is  that  modes  1-10  at  3515  km  range  contain  multiple 
arrivals,  rather  than  the  single,  dispersive  arrivals  that  are  associated  with  adia¬ 
batic  propagation.  This  has  important  implications  for  matched  field  processing 
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and  tomography  because  the  problem  of  associating  a  signal  with  a  particular  path 
through  the  ocean  becomes  difficult. 

Second,  many  of  the  arrivals  in  the  data  show  frequency-selective  fading  due  to  the 
destructive  interference  of  overlapping  multipaths.  There  is  evidence  of  frequency- 
coherent  arrivals  in  the  ATOC  receptions,  suggesting  that  the  STFT  processor  may 
be  resolving  some  individual  paths.  More  work  is  necessary  to  parameterize  the 
frequency-selectivity  of  the  channel  and  to  determine  how  often  coherent  arrivals 
can  be  expected. 

Third,  the  modal  multipath  structure  exhibits  significant  temporal  variability. 
Average  coherence  times4  for  the  modes  are  on  the  order  of  6-8  minutes,  which  is 
less  than  the  18  minute  source  transmission  time.  Coherence  appears  to  be  a  mild 
function  of  frequency,  but  not  a  function  of  mode  number  for  the  first  10  modes. 
Some  arrivals  are  temporally  stable  over  a  20-minute  reception,  but  are  not  usually 
detectable  in  the  next  reception,  four  hours  later.  The  results  of  this  chapter  indicate 
that  each  transmission  at  4-hour  intervals  essentially  samples  a  different  realization 
of  the  internal  wave  field. 

Fourth,  despite  the  complexity  of  the  mode  arrival  structure,  averaging  over 
source  transmissions  at  successive  4-hour  intervals  reveals  some  of  the  expected 
dispersion  characteristics  of  a  deep  water  channel.  An  analysis  of  three  statistics 
(leading  edge,  falling  edge,  and  centroid)  suggests  that  the  modes  retain  some  travel 
time  information  at  megameter  ranges.  The  results  of  this  chapter  indicate  that 
there  were  statistically  significant  changes  in  the  centroid  and  falling  edge  arrival 
times  over  the  5  months  of  data  analyzed  for  the  Hawaii  VLA.  On  average5,  the 
mode  arrival  times  decrease  by  0.4  seconds  between  December  and  March,  and  then 
increase  by  0.2  seconds  between  March  and  May.  The  RMS  errors  in  these  centroid 
measurements  are  20-30  ms  on  average,  on  the  same  order  as  fluctuations  of  ray 

i.e.,  the  time  at  which  the  correlation  coefficient  is  approximately  0.5 
5 Averaging  over  the  modes  in  the  75  Hz  bin 
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arrivals  and  pulse  termination  observed  by  other  researchers. 

Comparisons  of  the  ATOC  receptions  to  PE  simulations  through  internal  wave 
fields  at  1/2  Garrett-Munk  strength  support  these  conclusions.  Both  the  nature  of 
the  frequency-selective  fading  and  the  time  variations  of  the  modes  in  the  simulat¬ 
ed  data  qualitatively  agree  with  the  experimental  data.  Furthermore,  the  centroid 
statistics  demonstrate,  in  a  more  quantitative  way,  that  the  simulations  using  1/2 
GM  levels  are  in  agreement  with  the  data.  This  is  an  important  conclusion  because 
studies  of  the  ray-like  arrivals  indicate  internal  wave  strengths  on  the  order  of  1  /2 
GM. 

One  major  discrepancy  between  the  simulations  and  the  ATOC  data  is  that  the 
simulations  contain  a  sharp  axial  cutoff,  whereas  the  ATOC  data  does  not.  This  is 
likely  the  consequence  of  ignoring  the  bottom  interaction  near  the  Pioneer  Seamount 
source,  but  more  analysis  is  needed  to  verify  the  exact  coupling  mechanisms  that 
are  producing  the  complicated  excitation  pattern  in  the  low  modes.  In  addition, 
future  work  should  include  computing  the  mode  estimates  for  more  realizations  of 
the  internal  wave  field  so  that  temporal  statistics  for  the  simulations  can  be  compared 
to  the  measured  coherence  times  from  ATOC.  These  should  incorporate  the  effects 
of  the  sloping  bottom  near  the  source. 
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Chapter  5 


Conclusions  and  Future  Directions 


This  dissertation  has  developed  methods  of  estimating  the  modal  content  of  broad¬ 
band  receptions  and  used  those  methods  to  analyze  low-mode  arrivals  for  the  Acous¬ 
tic  Thermometry  of  Ocean  Climate  experiment.  This  chapter  will  summarize  the 
key  contributions  of  the  work  and  indicate  directions  for  future  research. 

The  first  contribution  of  this  thesis  is  a  general  framework  for  broadband  mode 
processing.  Most  previous  work  has  focused  on  narrowband  mode  estimation,  but 
recent  experiments  such  as  ATOC  require  wider  bandwidths.  This  thesis  has  de¬ 
veloped  a  mode  processor  based  on  the  short-time  Fourier  transform  and  addressed 
the  fundamental  issue  of  the  frequency  resolution  required  for  broadband  mode  es¬ 
timation.  A  key  advantage  of  the  STFT  approach  is  that  it  allows  study  of  the  the 
frequency  characteristics  of  individual  arrivals  (within  the  temporal  resolution  limits 
of  the  processor). 

Chapter  3  explored  the  time-  and  frequency-domain  characteristics  of  two  deter¬ 
ministic  modal  beamforming  algorithms:  the  matched  filter  and  the  pseudo-inverse 
filter.  A  detailed  performance  analysis  using  the  ATOC  array  configuration  showed 
that  it  is  possible  to  estimate  the  first  10  modes  using  a  pseudo-inverse  spatial  fil¬ 
ter,  and  a  frequency  bin  spacing  of  2.5  Hz.  The  latter  specification  implies  that  the 
temporal  resolution  of  the  processor  is  on  the  order  of  0.2  seconds  (depending  on  the 
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type  of  temporal  filter  used). 

The  second  major  contribution  of  this  work  is  an  investigation  of  the  low-mode 
arrivals  in  the  ATOC  Hawaii  data  set.  In  experiments  like  ATOC,  understanding 
the  axial  mode  arrivals  is  crucial  because  they  are  the  most  energetic  signals  at  long 
ranges.  The  analysis  presented  in  this  thesis  is  the  first  detailed  look  at  variability  of 
the  mode  signals  on  short  time  and  frequency  scales.  In  particular,  STFT  processing 
of  the  Hawaii  VLA  receptions  has  identified  several  important  features  of  the  mode 
signals  at  3515  km  range.  First,  each  low  mode  contains  series  of  arrivals,  rather 
than  the  single  dispersive  arrival  that  would  characterize  adiabatic  propagation.  This 
multipath  structure  exhibits  both  frequency-selective  fading  and  significant  temporal 
variability.  Average  coherence  times  are  on  the  order  of  6-8  minutes,  which  are 
somewhat  less  than  the  tens  of  minutes  predicted  by  Flatte  and  Stoughton  [72] 
and  the  12.7  minute  averaging  times  used  by  Worcester  et  al.  [74]  for  the  ATOC 
Engineering  Test  experiment.  (Note  that  since  some  peak  arrivals  in  the  STFT 
estimates  are  coherent  over  the  full  18.2  minute  transmission,  future  experiments 
should  include  longer  transmission  times  in  order  to  measure  the  actual  duration 
of  these  stable  peaks.)  'Given  this  high  degree  of  variability,  it  is  not  surprising 
that  receptions  at  4-hour  intervals  show  significant  differences.  Each  one  of  the 
transmissions  effectively  measures  a  different  realization  of  the  internal  wave  field. 
These  results  indicate  that  stochastic  methods  will  be  required  for  tomography  and 
matched  field  applications  that  use  the  mode  signals. 

The  mode  statistics  used  in  this  thesis  were  computed  by  averaging  over  groups 
of  receptions  at  4-hour  intervals.  The  leading  edge,  falling  edge,  and  centroid  of  the 
mode  arrivals  reveal  some  of  the  expected  dispersion  characteristics  of  a  deep  wa¬ 
ter  channel,  indicating  that  the  modes  retain  travel  time  information  at  megameter 
ranges.  Both  the  mode  centroids  and  the  falling  edges  showed  statistically  significant 
changes  in  mean  arrival  time  over  the  5  months  of  data  that  were  analyzed:  decreas¬ 
ing  arrival  times  from  December  to  March  followed  by  increasing  travel  times  from 
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March  to  May.  The  observed  trends  agree  with  results  for  the  ray  arrivals  over  the 
same  period  of  time.  As  noted  in  the  Science  article  by  the  ATOC  Consortium,  these 
trends  do  not  follow  the  seasonal  trends  in  the  surface  layer  and  are  likely  related  to 
a  subsurface  warming  near  the  Hawaii  array  [73].  The  RMS  errors  on  the  centroid 
and  falling  edge  statistics  are  on  the  order  of  20-30  ms  and  30-35  ms,  respectively. 

PE  simulations  of  the  propagation  through  internal  waves  at  1/2  Garrett-Munk 
strength  model  the  experimental  data  in  important  respects.  Since  1/2  GM  has 
proved  useful  in  explaining  some  of  the  internal  wave  effects  associated  with  the  ray 
arrivals,  this  is  a  useful  consistency  check.  The  primary  discrepancy  between  the 
simulations  and  the  ATOC  receptions  is  that  the  simulations  show  a  sharp  axial  cut¬ 
off,  whereas  the  experimental  data  does  not.  This  thesis  indicated  that  bathymetric 
coupling  near  the  Pioneer  Seamount  source  complicates  the  mode  excitation  pattern 
and  may  be  responsible  for  the  difference  between  simulated  and  experimental  data. 

In  terms  of  future  research,  this  thesis  has  laid  a  foundation  for  a  stochastic  chan¬ 
nel  model  that  is  useful  for  modal  tomography  and  other  applications.  In  pursuing 
this  ultimate  goal,  the  rich  set  of  data  provided  by  the  ATOC  experiment  merits 
further  study.  The  short-time  Fourier  framework  developed  here  could  be  extended 
to  include  methods  for  generating  a  time  series  from  the  time-varying  spectral  esti¬ 
mates  of  each  mode.  Fourier  synthesis  provides  a  useful  method  of  examining  which 
arrivals  add  coherently  across  frequency  bands.  Note  that  the  key  to  implementing 
the  synthesis  step  in  mode  processing  is  the  dispersion  correction.  As  indicated  in 
Chapter  3,  phase  shifts  across  frequency  due  to  dispersion  are  significant  at  mega¬ 
meter  ranges.  Since  scattering  may  produce  multiple  arrivals  in  a  single  mode,  each 
with  its  own  phase  characteristics,  the  correction  will  have  to  be  done  by  “steering” 
over  a  set  of  possible  corrections. 

Once  the  methods  for  synthesizing  modal  time  series  have  been  developed,  s- 
tatistics  of  the  broadband,  frequency-coherent  mode  arrivals  in  the  ATOC  data  can 
be  measured.  Time-of-arrival,  temporal  coherence,  and  intensity  fluctuations  are 
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especially  of  interest  since  comparable  statistics  exist  for  the  ray  arrivals.  The  re¬ 
lationship  of  rays  and  modes  propagating  through  random  internal  wave  fields  is 
currently  an  important  topic  of  research. 

In  addition  to  extending  the  STFT  framework  and  continuing  the  analysis  of 
the  Hawaii  data  set,  this  thesis  suggests  several  other  intriguing  topics  for  future 
research.  First,  the  problem  of  modal  beamforming  in  uncertain  environments  has 
been  raised,  but  not  solved,  by  this  thesis.  Based  on  the  degree  of  mismatch  between 
the  modes  for  the  environments  measured  at  deployment  and  recovery  of  the  ATOC 
arrays  (described  in  Chapter  2),  it  is  strongly  recommended  that  future  experiments 
include  regular  sampling  of  the  environment  at  the  array.  Measurements  of  the  local 
temperature  profile,  in  particular,  would  reduce  the  uncertainty  in  the  sound  speed 
profile,  thereby  reducing  the  potential  for  significant  modeshape  mismatch. 

Second,  an  in-depth  theoretical  analysis  of  the  modal  mismatch  problem  should 
be  undertaken.  Such  a  study  could  provide  useful  insights  about  the  robustness  of 
mode  filters  and  set  guidelines  for  the  required  temporal  and  spatial  sampling  of  the 
environment.  One  approach  to  this  problem  would  be  to  use  a  simple  dynamical 
model  for  the  mesoscale  sound  speed  fluctuations  and  relate  them  to  fluctuations 
in  the  modes  using  linear  perturbation  theory.  From  there  it  should  be  possible  to 
compute  performance  bounds  for  mode  estimators  under  a  variety  of  conditions  and 
to  explore  possible  strategies  for  mitigating  the  effects  of  environmental  uncertainty. 
This  is  a  relatively  open  problem,  although  there  is  a  large  body  of  related  work  in 
the  matched  field  processing  literature. 

Third,  a  study  of  the  excitation  of  the  axial  modes  by  a  bottom-mounted  source 
located  on  a  steep  slope  is  important.  The  obvious  starting  point  for  this  research 
is  to  model  the  propagation  down  the  slope  using  a  code,  such  as  range-dependent 
OASES,  that  accounts  for  the  elastic  properties  of  the  bottom.  In  addition  to  inves¬ 
tigating  the  forward  propagation  problem,  the  receptions  on  the  ATOC  VLA’s  from 
a  related  experiment  known  as  the  Alternate  Source  Test  (AST)  should  be  analyzed. 
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This  experiment  used  a  source  suspended  from  a  ship;  thus,  it  excited  the  low  modes 
directly.  A  comparison  of  the  AST  data  with  the  ATOC  data  might  provide  useful 
insights  into  the  effects  of  the  bathymetric  coupling  on  the  arrival  patterns.  Once 
the  excitation  pattern  of  the  low  modes  is  known,  an  interesting  question  to  ask  is 
whether  some  of  the  effects  of  the  bottom  interaction  can  be  removed  via  deconvo¬ 
lution  techniques.  Although  this  is  likely  to  be  a  difficult  problem,  a  solution  might 
provide  researchers  with  a  simpler  signal  for  use  in  studying  internal  wave  effects 
and/or  temperature  changes. 
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Appendix  A 


Acoustic  Mode  Perturbations 

As  discussed  in  Ch.  2,  underwater  acoustic  normal  modes  satisfy  a  second-order 
eigenvalue  equation: 

(f  ft2 

dz 2  ^  c2(z )  =  kmity&miz,  fi).  (A.l) 

Equation  A.l  assumes  unit  density  ( p(z )  =  1,  Vz),  which  is  reasonable  for  the  water¬ 
borne  modes.  These  modes  form  a  complete  set  of  orthogonal  basis  functions  and 
are  normalized  such  that 

f  Zmax 

Jo  <t>m(tt,z)(j)n(£l,z)dz  =  6(m  -  n),  (A.2) 

again  assuming  unit  density.  The  purpose  of  this  chapter  is  to  investigate  how  the 
mode  wavenumbers  and  shapes  are  affected  by  changes  in  using  linear  perturbation 
theory.  The  discussion  begins  with  a  brief  review  of  perturbation  theory.  Following 
that,  the  solution  to  the  problem  of  frequency  perturbations  is  straightforward. 
Linear  perturbation  theory  provides  a  method  of  solving  for  the  eigenvalues  and 
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eigenfunctions  of  a  perturbed  problem,  i.e., 

A 

— - — — - - s 

A  +  eAi  +  e2A2  H - ]  4>m  -  (A.3) 

in  terms  of  the  solutions  to  an  unperturbed  problem 

A<^ttl  l-l'mftm-  (A. 4) 

where  the  boundary  conditions  are  the  same  for  both  problems.  The  intent  of  this 
section  is  to  review  standard  results  found  in  a  number  of  references,  e.g.,  [75,  76], 
rather  than  to  provide  a  rigorous  derivation.  See  the  book  by  Rellich  for  a  thorough 
mathematical  treatment  of  the  problem  [77].  For  the  purposes  of  this  discussion,  the 
operator  A  has  distinct  eigenvalues  and  its  eigenvectors  form  a  complete  orthonormal 
(CON)  set.  Let  the  abbreviation  (x,  y )  represent  the  inner  product: 

{x,y)=  x(z)y(z)dz  (A.5) 

where  x  and  y  represent  arbitrary  functions  of  z.  The  following  results  rely  on  the 
standard  assumption  that  the  operator  A  is  self-adjoint,  i.e.,  ( Ax,y )  =  (x,  Ay). 

Assuming  that  A  may  be  written  in  terms  of  a  power  series  in  a  small  parameter 
e,  perturbation  theory  seeks  solutions  to  Eq.  A.3  of  the  form 

fim  =  ym  +  e/4^  +  e2^m  d - 

4>m  =  (f>m  +  +  e2^)  +  -  •  • 

Substituting  A. 6  into  Eq.  A.3  and  equating  terms  of  order  e  yields:1 

A(t>m  +  Al(j)m  =  +  y^(j)m  (A.7) 

xThe  e°  terms  represent  the  unperturbed  problem. 
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Taking  the  inner  product  of  both  sides  with  (f>m  and  solving  gives  the  first-order 
corrections  to  the  eigenvalues: 


9rn  — 


(A.8) 


Since  the  modes  of  the  unperturbed  problem  form  a  CON  set,  the  eigenfunction 
corrections  can  be  written  in  terms  of  the  modal  basis  set,  i.e., 

OO 

<f>m  =  E  9kl<l>k  (A-9) 

k=l 


where  g ^  is  the  coefficient  representing  the  contribution  of  the  kth.  mode  to  the  first- 
order  perturbation  of  the  mth  mode.  Substituting  Eq.  A. 9  into  Eq.  A. 3  and  taking 
the  inner  product  with  (f>k  ( k  ^  m)  gives  the  perturbation  expansion  coefficients: 


G(1)  =  < 
9fcm  \ 


(Ai4  4>k) 

9m  9k 


k  ^  m 


k  =  m 


(A-10) 


The  second-order  corrections  can  be  obtained  in  a  similar  manner. 


9rn  —  9km9mk(lJ'k  9m)  +  4m) 


(A.ll) 


and 


where 


«  =  E 


(A.12) 


fc=l 


9k 


9m  9k 


Y,9ki(9m  -  9k)  -  9m9kl  +  (A2 4m,  4k) 


L  k 


■5E(S^)2 


k  ^  m 


k  =  m 


■  (A. 13) 
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Consider  how  the  mode  parameters  vary  when  the  desired  frequency  is  equal  to 
w  =  to  +  A  to.  In  this  case,  the  perturbed  eigenvalue  problem  becomes 

d 2  u2  (  Au  /AuA2\1  ~  ~9~ 

d?  +  ^(z)  (1  +  2^r+  (v)  jj  =  (A.14) 

Using  the  following  definitions,  this  problem  fits  the  standard  form  discussed  above: 

Alj 
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